reverting r45444 to disentangle modules and make porting possible
[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 "AliHLTErrorGuard.h"
30 #include "TDirectory.h"
31 #include "TDatime.h"
32 #include "TString.h"
33 #include "TTree.h"
34 #include "TH1.h"
35 #include "TStopwatch.h"
36 #include "TUUID.h"
37 #include "TSystem.h"
38 #include "TRandom3.h"
39
40 /** ROOT macro for the implementation of ROOT specific class methods */
41 ClassImp(AliHLTTTreeProcessor)
42
43 AliHLTTTreeProcessor::AliHLTTTreeProcessor()
44                         : AliHLTProcessor(), 
45                           fDefinitions(),
46                           fTree(0),
47                           fMaxEntries(kMaxEntries),
48                           fPublishInterval(kInterval),
49                           fLastTime(0),
50                           fpEventTimer(NULL),
51                           fpCycleTimer(NULL),
52                           fMaxMemory(700000),
53                           fMaxEventTime(0),
54                           fNofEventsForce(0),
55                           fForcedEventsCount(0),
56                           fSkippedEventsCount(0),
57                           fNewEventsCount(0),
58                           fUniqueId(0),
59                           fIgnoreCycleTime(10),
60                           fCycleTimeFactor(1.0)
61 {
62   // see header file for class documentation
63   // or
64   // refer to README to build package
65   // or
66   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
67 }
68
69 const AliHLTUInt32_t AliHLTTTreeProcessor::fgkTimeScale=1000000; // ticks per second
70
71 AliHLTTTreeProcessor::~AliHLTTTreeProcessor()
72 {
73   // see header file for class documentation
74 }
75
76 AliHLTComponentDataType AliHLTTTreeProcessor::GetOutputDataType()
77 {
78   // get the component output data type
79   return kAliHLTDataTypeHistogram;
80 }
81
82 void AliHLTTTreeProcessor::GetOutputDataSize(unsigned long& constBase, double& inputMultiplier)
83 {
84   // get the output size estimator
85   //
86   if (!fDefinitions.size()) {
87     HLTError("Can not calculate output data size, no histogram definitions were provided");
88     return;
89   }
90
91   constBase = 0;
92   for (list_const_iterator i = fDefinitions.begin(); i != fDefinitions.end(); ++i)
93     constBase += i->GetSize();
94
95   inputMultiplier = 1.;
96 }
97
98 int AliHLTTTreeProcessor::DoInit(int argc, const char** argv)
99 {
100   // init component
101   // ask child to create the tree.
102   int iResult = 0;
103
104   // calculating a unique id from the hostname and process id
105   // used for identifying output of multiple components
106   TUUID guid = GenerateGUID();
107   union
108   {
109     UChar_t buf[16];
110     UInt_t bufAsInt[4];
111   };
112   guid.GetUUID(buf);
113   fUniqueId = bufAsInt[0];
114   
115
116   if (!fTree) {
117     std::auto_ptr<TTree> ptr(CreateTree(argc, argv));
118     if (ptr.get()) {
119       //Stage 1: default initialization.
120       ptr->SetDirectory(0);
121       //"Default" (for derived component) histograms.
122       FillHistogramDefinitions();
123       //Default values.
124       fMaxEntries = kMaxEntries;
125       fPublishInterval = kInterval;
126       fLastTime = 0;
127       //Stage 2: OCDB.
128       TString cdbPath("HLT/ConfigHLT/");
129       cdbPath += GetComponentID();
130       //
131       iResult = ConfigureFromCDBTObjString(cdbPath);
132       //
133       if (iResult < 0)
134         return iResult;
135       //Stage 3: command line arguments.
136       if (argc && (iResult = ConfigureFromArgumentString(argc, argv)) < 0)
137         return iResult;
138
139       ptr->SetCircular(fMaxEntries);
140       fTree = ptr.release();
141     } else //No way to process error correctly - error is unknown here.
142       return -EINVAL;
143   } else {
144     HLTError("fTree pointer must be null before DoInit call");
145     return -EINVAL;
146   }
147
148   if (iResult>=0 && fMaxEventTime>0) {
149     fpEventTimer=new TStopwatch;
150     if (fpEventTimer) {
151       fpEventTimer->Reset();
152     }
153     fpCycleTimer=new TStopwatch;
154     if (fpCycleTimer) {
155       fpCycleTimer->Reset();
156     }
157   }
158   fSkippedEventsCount=0;
159
160   return iResult;
161 }
162
163 int AliHLTTTreeProcessor::DoDeinit()
164 {
165   // cleanup component
166   delete fTree;
167   fTree = 0;
168   fDefinitions.clear();
169
170   if (fpEventTimer) delete fpEventTimer;
171   fpEventTimer=NULL;
172   if (fpCycleTimer) delete fpCycleTimer;
173   fpCycleTimer=NULL;
174
175   return 0;
176 }
177
178 int AliHLTTTreeProcessor::DoEvent(const AliHLTComponentEventData& evtData, AliHLTComponentTriggerData& trigData)
179 {
180   //Process event and publish histograms.
181   AliHLTUInt32_t eventType=0;
182   if (!IsDataEvent(&eventType) && eventType!=gkAliEventTypeEndOfRun) return 0;
183
184   //I'm pretty sure, that if fTree == 0 (DoInit failed) DoEvent is not called.
185   //But interface itself does not force you to call DoInit before DoEvent, so,
186   //I make this check explicit.
187   if (!fTree) {
188     HLTError("fTree is a null pointer, try to call AliHLTTTreeProcessor::DoInit first.");
189     return -EINVAL;//-ENULLTREE? :)
190   }
191
192   AliHLTUInt32_t averageEventTime=0;
193   AliHLTUInt32_t averageCycleTime=0;
194
195   int fillingtime=0;
196   int publishtime=0;
197   bool bDoFilling=true;
198   bool bDoPublishing=false;
199   const int cycleResetInterval=1000;
200   if (fpEventTimer && fpCycleTimer) {
201     averageEventTime=(fpEventTimer->RealTime()*fgkTimeScale)/(GetEventCount()+1);
202     fillingtime=fpEventTimer->RealTime()*fgkTimeScale;
203     publishtime=fillingtime;
204     fpEventTimer->Start(kFALSE);
205     fpCycleTimer->Stop();
206     averageCycleTime=(fpCycleTimer->RealTime()*fgkTimeScale)/((GetEventCount()%cycleResetInterval)+1);
207     // adapt processing to 3/4 of the max time
208     bDoFilling=4*averageEventTime<3*fMaxEventTime ||
209       (averageEventTime<fCycleTimeFactor*averageCycleTime && fpCycleTimer->RealTime()>fIgnoreCycleTime);
210     if (fNofEventsForce>0 && fForcedEventsCount<fNofEventsForce) {
211       fForcedEventsCount++;
212       bDoFilling=true;
213     }
214   }
215
216   // FIXME: there is still an unclear increase in memory consumption, even if the number of entries
217   // in the tree is restricted. Valgrind studies did not show an obvious memory leak. This is likely
218   // to be caused by something deep in the Root TTree functionality and needs to be studied in detail.
219   ProcInfo_t ProcInfo;
220   gSystem->GetProcInfo(&ProcInfo);
221   if (ProcInfo.fMemResident>fMaxMemory) bDoFilling=false;
222
223   // process input data blocks and fill the tree
224   int iResult = 0;
225   if (eventType!=gkAliEventTypeEndOfRun) {
226     if (bDoFilling) {iResult=FillTree(fTree, evtData, trigData); fNewEventsCount++;}
227     else fSkippedEventsCount++;
228   }
229   if (fpEventTimer) {
230     fpEventTimer->Stop();
231     fillingtime=fpEventTimer->RealTime()*fgkTimeScale-fillingtime;
232     if (fillingtime<0) fillingtime=0;
233     fpEventTimer->Start(kFALSE);
234   }
235
236   if (iResult < 0) {
237     ALIHLTERRORGUARD(5, "FillTree failed with %d, first event %d", iResult, GetEventCount());
238     return iResult;
239   }
240
241   const TDatime time;
242
243   if (( time.Get() - fLastTime > fPublishInterval && fNewEventsCount>0) ||
244       eventType==gkAliEventTypeEndOfRun) {
245     if ((bDoPublishing=fLastTime>0)) { // publish earliest after the first interval but set the timer
246
247     for (list_const_iterator i = fDefinitions.begin(); i != fDefinitions.end(); ++i) {
248       if (TH1* h = CreateHistogram(*i)) {
249         //I do not care about errors here - since I'm not able
250         //to rollback changes.
251         // TODO: in case of -ENOSPC et the size of the last object by calling
252         // GetLastObjectSize() and accumulate the necessary output buffer size
253         PushBack(h, GetOriginDataType(), GetDataSpec());
254         delete h;
255       }
256     }
257     unsigned eventcount=GetEventCount()+1;
258     HLTBenchmark("publishing %d histograms, %d entries in tree, %d new events since last publishing, accumulated %d of %d events (%.1f%%)", fDefinitions.size(), fTree->GetEntriesFast(), fNewEventsCount, eventcount-fSkippedEventsCount, eventcount, eventcount>0?(100*float(eventcount-fSkippedEventsCount)/eventcount):0);
259     fNewEventsCount=0;
260     HLTBenchmark("current memory usage %d %d", ProcInfo.fMemResident, ProcInfo.fMemVirtual);
261     }
262
263     fLastTime=time.Get();
264     if (fLastTime==0) {
265       // choose a random offset at beginning to equalize traffic for multiple instances
266       // of the component
267       gRandom->SetSeed(fUniqueId);
268       fLastTime-=gRandom->Integer(fPublishInterval);
269     }
270   }
271
272   if (fpEventTimer) {
273     fpEventTimer->Stop();
274     publishtime=fpEventTimer->RealTime()*fgkTimeScale-publishtime;
275     if (publishtime>fillingtime) publishtime-=fillingtime;
276     else publishtime=0;
277
278     averageEventTime=(fpEventTimer->RealTime()*fgkTimeScale)/(GetEventCount()+1);
279
280     // info output once every 5 seconds
281     static UInt_t lastTime=0;
282     if (time.Get()-lastTime>5 ||
283         eventType==gkAliEventTypeEndOfRun ||
284         bDoPublishing) {
285       lastTime=time.Get();
286       unsigned eventcount=GetEventCount()+1;
287       HLTBenchmark("filling time %d us, publishing time %d, average total processing time %d us, cycle time %d us, accumulated %d of %d events (%.1f%%)", fillingtime, publishtime, averageEventTime, averageCycleTime, eventcount-fSkippedEventsCount, eventcount, eventcount>0?(100*float(eventcount-fSkippedEventsCount)/eventcount):0);
288     }
289   }
290   if (fpCycleTimer) {
291     bool bReset=(GetEventCount()%cycleResetInterval)==0;
292     fpCycleTimer->Start(bReset);
293   }
294
295   return iResult;
296 }
297
298 int AliHLTTTreeProcessor::ScanConfigurationArgument(int argc, const char** argv)
299 {
300   // scan one argument and its parameters from the list
301   // return number of processed entries.
302   // possible arguments: 
303   // -maxentries number
304   // -interval number
305   // -histogram name -size number -expression expression [-cut expression ][-opt option]
306   // As soon as "-histogram" found, -size and -expression and -outtype are required, 
307   // cut and option can be omitted.
308   if (argc <= 0)
309     return 0;
310
311   std::list<AliHLTHistogramDefinition> newDefs;
312   AliHLTHistogramDefinition def;
313
314   int i = 0;
315   int maxEntries = 0;
316
317   while (i < argc) {
318     const TString argument(argv[i]);
319
320     if (argument.CompareTo("-maxentries") == 0) { //1. Max entries argument for TTree.
321       if (i + 1 == argc) {
322         HLTError("Numeric value for '-maxentries' is expected");
323         return -EPROTO;
324       }
325       //Next must be a number.
326       //TString returns 0 (number) even if string contains non-numeric symbols.
327       maxEntries = TString(argv[i + 1]).Atoi();
328       if (maxEntries <= 0) {
329         HLTError("Bad value for '-maxentries': %d", maxEntries);
330         return -EPROTO;
331       }
332   
333       i += 2;
334     } else if (argument.CompareTo("-interval") == 0) { //2. Interval argument for publishing.
335       if (i + 1 == argc) {
336         HLTError("Numeric value for '-interval' is expected");
337         return -EPROTO;
338       }
339
340       const Int_t interval = TString(argv[i + 1]).Atoi();
341       if (interval < 0) {
342         HLTError("Bad value for '-interval' argument: %d", interval);
343         return -EPROTO;
344       }
345
346       fPublishInterval = interval;
347
348       i += 2;
349     } else if (argument.CompareTo("-maxeventtime") == 0) { // max average processing time in us
350       if (i + 1 == argc) {
351         HLTError("Numeric value for '-maxeventtime' is expected");
352         return -EPROTO;
353       }
354
355       const Int_t time = TString(argv[i + 1]).Atoi();
356       if (time < 0) {
357         HLTError("Bad value for '-maxeventtime' argument: %d", time);
358         return -EPROTO;
359       }
360
361       fMaxEventTime = time;
362
363       i += 2;
364     } else if (argument.CompareTo("-forced-events") == 0) { // number of forced events
365       if (i + 1 == argc) {
366         HLTError("Numeric value for '-forced-events' is expected");
367         return -EPROTO;
368       }
369
370       const Int_t count = TString(argv[i + 1]).Atoi();
371       if (count < 0) {
372         HLTError("Bad value for '-forced-events' argument: %d", count);
373         return -EPROTO;
374       }
375
376       fNofEventsForce = count;
377       fForcedEventsCount=0;
378
379       i += 2;
380     } else if (argument.CompareTo("-ignore-cycletime") == 0) { // ignore cycle time for n sec
381       if (i + 1 == argc) {
382         HLTError("Numeric value for '-ignore-cycletime' is expected");
383         return -EPROTO;
384       }
385
386       const Int_t time = TString(argv[i + 1]).Atoi();
387       if (time < 0) {
388         HLTError("Bad value for '-ignore-cycletime' argument: %d", time);
389         return -EPROTO;
390       }
391
392       fIgnoreCycleTime = time;
393       i += 2;
394     } else if (argument.CompareTo("-maxmemory") == 0) { // maximum of memory in kByte to be used by the component
395       if (i + 1 == argc) {
396         HLTError("Numeric value for '-maxmemory' is expected");
397         return -EPROTO;
398       }
399
400       const Int_t mem = TString(argv[i + 1]).Atoi();
401       if (mem < 0) {
402         HLTError("Bad value for '-maxmemory' argument: %d", time);
403         return -EPROTO;
404       }
405
406       fMaxMemory = mem;
407       i += 2;
408     } else if (argument.CompareTo("-cycletime-factor") == 0) { // weight factor for cycle time
409       if (i + 1 == argc) {
410         HLTError("Numeric value for '-cycletime-factor' is expected");
411         return -EPROTO;
412       }
413
414       const Float_t factor = TString(argv[i + 1]).Atof();
415       if (factor < 0) {
416         HLTError("Bad value for '-cycletime-factor' argument: %f", factor);
417         return -EPROTO;
418       }
419
420       fCycleTimeFactor = factor;
421       i += 2;
422     } else if (argument.CompareTo("-histogram") == 0) { //3. Histogramm definition.
423       const int nParsed = ParseHistogramDefinition(argc, argv, i, def);
424       if (!nParsed)
425         return -EPROTO;
426
427       newDefs.push_back(def);
428
429       i += nParsed;   
430     } else {
431       HLTError("Unknown argument %s", argument.Data());
432       return -EPROTO;
433     }
434   }
435
436   if (maxEntries != fMaxEntries) {
437     fMaxEntries = maxEntries;
438     if (fTree) {
439       fTree->Reset();
440       fTree->SetCircular(fMaxEntries);
441     }
442   }
443
444   if (newDefs.size())
445     fDefinitions.swap(newDefs);
446
447   return i;
448 }
449
450 TH1* AliHLTTTreeProcessor::CreateHistogram(const AliHLTHistogramDefinition& d)
451 {
452
453   // create a histogram from the tree
454   if (!fTree) {
455     HLTError("fTree is a null pointer, try to call AliHLTTTreeProcessor::DoInit first.");
456     return 0;
457   }
458
459   TString histName(d.GetName());
460   if (!histName.Contains("(")) {
461     //Without number of bins, the histogram will be "fixed"
462     //and most of values can go to underflow/overflow bins,
463     //since kCanRebin will be false.
464     histName += TString::Format("(%d)", Int_t(kDefaultNBins));
465   }
466
467   const Long64_t rez = fTree->Project(histName.Data(), d.GetExpression().Data(), d.GetCut().Data(), d.GetDrawOption().Data());
468
469   if (rez == -1) {
470     HLTError("TTree::Project failed");
471     return 0;
472   }
473
474   //Now, cut off the binning part of a name
475   histName = histName(0, histName.Index("("));
476   TH1 * hist = dynamic_cast<TH1*>(gDirectory->Get(histName.Data()));
477   if (!hist) {
478     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()));
479     HLTError(msg.Data());
480   }else if (d.GetDrawOption().Length()) {
481     hist->SetOption(d.GetDrawOption().Data());
482   }
483
484   return hist;
485 }
486
487 int AliHLTTTreeProcessor::ParseHistogramDefinition(int argc, const char** argv, int pos, AliHLTHistogramDefinition& dst)const
488 {
489   //Histogram-definition:
490   //    -histogram name -size number -expression expression [-cut expression][-opt option]
491
492   //at pos we have '-histogram', at pos + 1 must be the name.
493   if (pos + 1 == argc) {
494     HLTError("Bad histogram definition, histogram name is expected");
495     return 0;
496   }
497
498   dst.SetName(argv[pos + 1]);
499   pos += 2;
500   
501   //At pos must be '-size', and number at pos + 1.
502   if (pos == argc || TString(argv[pos]).CompareTo("-size")) {
503     HLTError("Bad histogram definition, '-size' is expected");
504     return 0;
505   }
506
507   if (pos + 1 == argc) {
508     HLTError("Bad histogram definition, size is expected");
509     return 0;
510   }
511
512   dst.SetSize(TString(argv[pos + 1]).Atoi());
513   if (dst.GetSize() <= 0) {
514     HLTError("Bad histogram definition, positive size is required");
515     return 0;
516   }
517
518   pos += 2;
519   //At pos must be '-expression', and expression at pos + 1. 
520   if (pos == argc || TString(argv[pos]).CompareTo("-expression")) {
521     HLTError("Bad histogram definition, '-expression' is expected");
522     return 0;
523   }
524
525   if (pos + 1 == argc) {
526     HLTError("Bad histogram definition, expression is expected");
527     return 0;
528   }
529
530   dst.SetExpression(argv[pos + 1]);
531   pos += 2;
532
533   int processed = 6;
534   dst.SetCut("");
535   dst.SetDrawOption("");
536
537   //remaining options can be the cut and Draw option.
538   //cut must be first.
539   if (pos + 1 >= argc)
540     return processed;
541
542   if (TString(argv[pos]).CompareTo("-cut") == 0) {
543     dst.SetCut(argv[pos + 1]);
544     pos += 2;
545     processed += 2;
546   }
547
548   if (pos + 1 >= argc)
549     return processed;
550
551   if (TString(argv[pos]).CompareTo("-opt") == 0) {
552     dst.SetDrawOption(argv[pos + 1]);
553     processed += 2;
554   }
555
556   return processed;
557 }