2c482605c004a543db69a21546504d1d01083f75
[u/mrichter/AliRoot.git] / HLT / BASE / AliHLTTTreeProcessor.cxx
1 // $Id$
2
3 //**************************************************************************
4 //* This file is property of and copyright by the ALICE HLT Project        * 
5 //* ALICE Experiment at CERN, All rights reserved.                         *
6 //*                                                                        *
7 //* Primary Authors: Timur Pocheptsov <Timur.Pocheptsov@cern.ch>           *
8 //*                  Matthias Richter <Matthias.Richter@cern.ch>
9 //*                  for The ALICE HLT Project.                            *
10 //*                                                                        *
11 //* Permission to use, copy, modify and distribute this software and its   *
12 //* documentation strictly for non-commercial purposes is hereby granted   *
13 //* without fee, provided that the above copyright notice appears in all   *
14 //* copies and that both the copyright notice and this permission notice   *
15 //* appear in the supporting documentation. The authors make no claims     *
16 //* about the suitability of this software for any purpose. It is          *
17 //* provided "as is" without express or implied warranty.                  *
18 //**************************************************************************
19
20 /// @file   AliHLTTTreeProcessor.cxx
21 /// @author Timur Pocheptsov, Matthias Richter
22 /// @date   05.07.2010
23 /// @brief  Generic component for data collection in a TTree
24
25 #include <cerrno>
26 #include <memory>
27
28 #include "AliHLTTTreeProcessor.h"
29 #include "TDirectory.h"
30 #include "TDatime.h"
31 #include "TString.h"
32 #include "TTree.h"
33 #include "TH1.h"
34 #include "TStopwatch.h"
35
36 /** ROOT macro for the implementation of ROOT specific class methods */
37 ClassImp(AliHLTTTreeProcessor)
38
39 AliHLTTTreeProcessor::AliHLTTTreeProcessor()
40                         : AliHLTProcessor(), 
41                           fDefinitions(),
42                           fTree(0),
43                           fMaxEntries(kMaxEntries),
44                           fPublishInterval(kInterval),
45                           fLastTime(0),
46                           fpEventTimer(NULL),
47                           fpCycleTimer(NULL),
48                           fMaxEventTime(0),
49                           fNofEventsForce(0),
50                           fForcedEventsCount(0),
51                           fSkippedEventsCount(0)
52 {
53   // see header file for class documentation
54   // or
55   // refer to README to build package
56   // or
57   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
58 }
59
60 const AliHLTUInt32_t AliHLTTTreeProcessor::fgkTimeScale=1000000; // ticks per second
61
62 AliHLTTTreeProcessor::~AliHLTTTreeProcessor()
63 {
64   // see header file for class documentation
65 }
66
67 AliHLTComponentDataType AliHLTTTreeProcessor::GetOutputDataType()
68 {
69   // get the component output data type
70   return kAliHLTDataTypeHistogram;
71 }
72
73 void AliHLTTTreeProcessor::GetOutputDataSize(unsigned long& constBase, double& inputMultiplier)
74 {
75   // get the output size estimator
76   //
77   if (!fDefinitions.size()) {
78     HLTError("Can not calculate output data size, no histogram definitions were provided");
79     return;
80   }
81
82   constBase = 0;
83   for (list_const_iterator i = fDefinitions.begin(); i != fDefinitions.end(); ++i)
84     constBase += i->GetSize();
85
86   inputMultiplier = 1.;
87 }
88
89 int AliHLTTTreeProcessor::DoInit(int argc, const char** argv)
90 {
91   // init component
92   // ask child to create the tree.
93   int iResult = 0;
94
95   if (!fTree) {
96     std::auto_ptr<TTree> ptr(CreateTree(argc, argv));
97     if (ptr.get()) {
98       //Stage 1: default initialization.
99       ptr->SetDirectory(0);
100       //"Default" (for derived component) histograms.
101       FillHistogramDefinitions();
102       //Default values.
103       fMaxEntries = kMaxEntries;
104       fPublishInterval = kInterval;
105       fLastTime = 0;
106       //Stage 2: OCDB.
107       TString cdbPath("HLT/ConfigHLT/");
108       cdbPath += GetComponentID();
109       //
110       iResult = ConfigureFromCDBTObjString(cdbPath);
111       //
112       if (iResult < 0)
113         return iResult;
114       //Stage 3: command line arguments.
115       if (argc && (iResult = ConfigureFromArgumentString(argc, argv)) < 0)
116         return iResult;
117
118       ptr->SetCircular(fMaxEntries);
119       fTree = ptr.release();
120     } else //No way to process error correctly - error is unknown here.
121       return -EINVAL;
122   } else {
123     HLTError("fTree pointer must be null before DoInit call");
124     return -EINVAL;
125   }
126
127   if (iResult>=0 && fMaxEventTime>0) {
128     fpEventTimer=new TStopwatch;
129     if (fpEventTimer) {
130       fpEventTimer->Reset();
131     }
132     fpCycleTimer=new TStopwatch;
133     if (fpCycleTimer) {
134       fpCycleTimer->Reset();
135     }
136   }
137   fSkippedEventsCount=0;
138
139   return iResult;
140 }
141
142 int AliHLTTTreeProcessor::DoDeinit()
143 {
144   // cleanup component
145   delete fTree;
146   fTree = 0;
147   fDefinitions.clear();
148
149   if (fpEventTimer) delete fpEventTimer;
150   fpEventTimer=NULL;
151   if (fpCycleTimer) delete fpCycleTimer;
152   fpCycleTimer=NULL;
153
154   return 0;
155 }
156
157 int AliHLTTTreeProcessor::DoEvent(const AliHLTComponentEventData& evtData, AliHLTComponentTriggerData& trigData)
158 {
159   //Process event and publish histograms.
160   AliHLTUInt32_t eventType=0;
161   if (!IsDataEvent(&eventType) && eventType!=gkAliEventTypeEndOfRun) return 0;
162
163   //I'm pretty sure, that if fTree == 0 (DoInit failed) DoEvent is not called.
164   //But interface itself does not force you to call DoInit before DoEvent, so,
165   //I make this check explicit.
166   if (!fTree) {
167     HLTError("fTree is a null pointer, try to call AliHLTTTreeProcessor::DoInit first.");
168     return -EINVAL;//-ENULLTREE? :)
169   }
170
171   AliHLTUInt32_t averageEventTime=0;
172   AliHLTUInt32_t averageCycleTime=0;
173
174   AliHLTUInt32_t proctime=0;
175   bool bDoFilling=true;
176   if (fpEventTimer && fpCycleTimer) {
177     averageEventTime=(fpEventTimer->RealTime()*fgkTimeScale)/(GetEventCount()+1);
178     proctime=fpEventTimer->RealTime()*fgkTimeScale;
179     fpEventTimer->Start(kFALSE);
180     fpCycleTimer->Stop();
181     averageCycleTime=(fpCycleTimer->RealTime()*fgkTimeScale)/(GetEventCount()+1);
182     // adapt processing to 3/4 of the max time
183     bDoFilling=4*averageEventTime<3*fMaxEventTime || averageEventTime<averageCycleTime;
184     if (fNofEventsForce>0 && fForcedEventsCount<fNofEventsForce) {
185       fForcedEventsCount++;
186       bDoFilling=true;
187     }
188   }
189
190   // process input data blocks and fill the tree
191   int iResult = 0;
192   if (eventType!=gkAliEventTypeEndOfRun) {
193     if (bDoFilling) iResult=FillTree(fTree, evtData, trigData);
194     else fSkippedEventsCount++;
195   }
196
197   if (iResult < 0)
198     return iResult;
199
200   const TDatime time;
201
202   if ( time.Get() - fLastTime > fPublishInterval ||
203       eventType==gkAliEventTypeEndOfRun) {
204     for (list_const_iterator i = fDefinitions.begin(); i != fDefinitions.end(); ++i) {
205       if (TH1* h = CreateHistogram(*i)) {
206         //I do not care about errors here - since I'm not able
207         //to rollback changes.
208         // TODO: in case of -ENOSPC et the size of the last object by calling
209         // GetLastObjectSize() and accumulate the necessary output buffer size
210         PushBack(h, GetOriginDataType(), GetDataSpec());
211       }
212     }
213
214     fLastTime = time.Get();
215   }
216
217   if (fpEventTimer) {
218     fpEventTimer->Stop();
219     proctime=fpEventTimer->RealTime()*fgkTimeScale-proctime;
220     averageEventTime=(fpEventTimer->RealTime()*fgkTimeScale)/(GetEventCount()+1);
221
222     // info output once every 5 seconds
223     static UInt_t lastTime=0;
224     if (time.Get()-lastTime>5 ||
225       eventType==gkAliEventTypeEndOfRun) {
226       lastTime=time.Get();
227       unsigned eventcount=GetEventCount();
228       HLTBenchmark("event time %d us, average time %d us, cycle time %d us, accumulated %d of %d events (%.1f%%)", proctime, averageEventTime, averageCycleTime, eventcount-fSkippedEventsCount, eventcount, eventcount>0?(100*float(eventcount-fSkippedEventsCount)/eventcount):0);
229     }
230   }
231   if (fpCycleTimer) {
232     fpCycleTimer->Start(kFALSE);
233   }
234
235   return iResult;
236 }
237
238 int AliHLTTTreeProcessor::ScanConfigurationArgument(int argc, const char** argv)
239 {
240   // scan one argument and its parameters from the list
241   // return number of processed entries.
242   // possible arguments: 
243   // -maxentries number
244   // -interval number
245   // -histogram name -size number -expression expression [-cut expression ][-opt option]
246   // As soon as "-histogram" found, -size and -expression and -outtype are required, 
247   // cut and option can be omitted.
248   if (argc <= 0)
249     return 0;
250
251   std::list<AliHLTHistogramDefinition> newDefs;
252   AliHLTHistogramDefinition def;
253
254   int i = 0;
255   int maxEntries = 0;
256
257   while (i < argc) {
258     const TString argument(argv[i]);
259
260     if (argument.CompareTo("-maxentries") == 0) { //1. Max entries argument for TTree.
261       if (i + 1 == argc) {
262         HLTError("Numeric value for '-maxentries' is expected");
263         return -EPROTO;
264       }
265       //Next must be a number.
266       //TString returns 0 (number) even if string contains non-numeric symbols.
267       maxEntries = TString(argv[i + 1]).Atoi();
268       if (maxEntries <= 0) {
269         HLTError("Bad value for '-maxentries': %d", maxEntries);
270         return -EPROTO;
271       }
272   
273       i += 2;
274     } else if (argument.CompareTo("-interval") == 0) { //2. Interval argument for publishing.
275       if (i + 1 == argc) {
276         HLTError("Numeric value for '-interval' is expected");
277         return -EPROTO;
278       }
279
280       const Int_t interval = TString(argv[i + 1]).Atoi();
281       if (interval <= 0) {
282         HLTError("Bad value for '-interval' argument: %d", interval);
283         return -EPROTO;
284       }
285
286       fPublishInterval = interval;
287
288       i += 2;
289     } else if (argument.CompareTo("-maxeventtime") == 0) { // max average processing time in us
290       if (i + 1 == argc) {
291         HLTError("Numeric value for '-maxeventtime' is expected");
292         return -EPROTO;
293       }
294
295       const Int_t time = TString(argv[i + 1]).Atoi();
296       if (time <= 0) {
297         HLTError("Bad value for '-maxeventtime' argument: %d", time);
298         return -EPROTO;
299       }
300
301       fMaxEventTime = time;
302
303       i += 2;
304     } else if (argument.CompareTo("-forced-events") == 0) { // number of forced events
305       if (i + 1 == argc) {
306         HLTError("Numeric value for '-forced-events' is expected");
307         return -EPROTO;
308       }
309
310       const Int_t count = TString(argv[i + 1]).Atoi();
311       if (count <= 0) {
312         HLTError("Bad value for '-forced-events' argument: %d", count);
313         return -EPROTO;
314       }
315
316       fNofEventsForce = count;
317       fForcedEventsCount=0;
318
319       i += 2;
320     } else if (argument.CompareTo("-histogram") == 0) { //3. Histogramm definition.
321       const int nParsed = ParseHistogramDefinition(argc, argv, i, def);
322       if (!nParsed)
323         return -EPROTO;
324
325       newDefs.push_back(def);
326
327       i += nParsed;   
328     } else {
329       HLTError("Unknown argument %s", argument.Data());
330       return -EPROTO;
331     }
332   }
333
334   if (maxEntries != fMaxEntries) {
335     fMaxEntries = maxEntries;
336     if (fTree) {
337       fTree->Reset();
338       fTree->SetCircular(fMaxEntries);
339     }
340   }
341
342   if (newDefs.size())
343     fDefinitions.swap(newDefs);
344
345   return i;
346 }
347
348 TH1* AliHLTTTreeProcessor::CreateHistogram(const AliHLTHistogramDefinition& d)
349 {
350
351   // create a histogram from the tree
352   if (!fTree) {
353     HLTError("fTree is a null pointer, try to call AliHLTTTreeProcessor::DoInit first.");
354     return 0;
355   }
356
357   TString histName(d.GetName());
358   if (!histName.Contains("(")) {
359     //Without number of bins, the histogram will be "fixed"
360     //and most of values can go to underflow/overflow bins,
361     //since kCanRebin will be false.
362     histName += TString::Format("(%d)", Int_t(kDefaultNBins));
363   }
364
365   const Long64_t rez = fTree->Project(histName.Data(), d.GetExpression().Data(), d.GetCut().Data(), d.GetDrawOption().Data());
366
367   if (rez == -1) {
368     HLTError("TTree::Project failed");
369     return 0;
370   }
371
372   //Now, cut off the binning part of a name
373   histName = histName(0, histName.Index("("));
374   TH1 * hist = dynamic_cast<TH1*>(gDirectory->Get(histName.Data()));
375   if (!hist) {
376     const TString msg(Form("Hist %s is a null pointer, selection was %s, strange name or hist's type\n", histName.Data(), d.GetExpression().Data()));
377     HLTError(msg.Data());
378   }else if (d.GetDrawOption().Length()) {
379     hist->SetOption(d.GetDrawOption().Data());
380   }
381
382   return hist;
383 }
384
385 int AliHLTTTreeProcessor::ParseHistogramDefinition(int argc, const char** argv, int pos, AliHLTHistogramDefinition& dst)const
386 {
387   //Histogram-definition:
388   //    -histogram name -size number -expression expression [-cut expression][-opt option]
389
390   //at pos we have '-histogram', at pos + 1 must be the name.
391   if (pos + 1 == argc) {
392     HLTError("Bad histogram definition, histogram name is expected");
393     return 0;
394   }
395
396   dst.SetName(argv[pos + 1]);
397   pos += 2;
398   
399   //At pos must be '-size', and number at pos + 1.
400   if (pos == argc || TString(argv[pos]).CompareTo("-size")) {
401     HLTError("Bad histogram definition, '-size' is expected");
402     return 0;
403   }
404
405   if (pos + 1 == argc) {
406     HLTError("Bad histogram definition, size is expected");
407     return 0;
408   }
409
410   dst.SetSize(TString(argv[pos + 1]).Atoi());
411   if (dst.GetSize() <= 0) {
412     HLTError("Bad histogram definition, positive size is required");
413     return 0;
414   }
415
416   pos += 2;
417   //At pos must be '-expression', and expression at pos + 1. 
418   if (pos == argc || TString(argv[pos]).CompareTo("-expression")) {
419     HLTError("Bad histogram definition, '-expression' is expected");
420     return 0;
421   }
422
423   if (pos + 1 == argc) {
424     HLTError("Bad histogram definition, expression is expected");
425     return 0;
426   }
427
428   dst.SetExpression(argv[pos + 1]);
429   pos += 2;
430
431   int processed = 6;
432   dst.SetCut("");
433   dst.SetDrawOption("");
434
435   //remaining options can be the cut and Draw option.
436   //cut must be first.
437   if (pos + 1 >= argc)
438     return processed;
439
440   if (TString(argv[pos]).CompareTo("-cut") == 0) {
441     dst.SetCut(argv[pos + 1]);
442     pos += 2;
443     processed += 2;
444   }
445
446   if (pos + 1 >= argc)
447     return processed;
448
449   if (TString(argv[pos]).CompareTo("-opt") == 0) {
450     dst.SetDrawOption(argv[pos + 1]);
451     processed += 2;
452   }
453
454   return processed;
455 }