FairRoot
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
FairMonitor.cxx
Go to the documentation of this file.
1 /********************************************************************************
2  * Copyright (C) 2014 GSI Helmholtzzentrum fuer Schwerionenforschung GmbH *
3  * *
4  * This software is distributed under the terms of the *
5  * GNU Lesser General Public Licence (LGPL) version 3, *
6  * copied verbatim in the file "LICENSE" *
7  ********************************************************************************/
8 /*
9  * FairMonitor.cxx
10  *
11  * Created on: Aug 01, 2015
12  * Author: r.karabowicz
13  */
14 
15 #include "FairMonitor.h"
16 
17 #include "FairLogger.h"
18 #include "FairSystemInfo.h"
19 
20 #include <TArrow.h>
21 #include <TAxis.h>
22 #include <TBox.h>
23 #include <TCanvas.h>
24 #include <TCollection.h>
25 #include <TFile.h>
26 #include <TH1.h>
27 #include <TLatex.h>
28 #include <TLegend.h>
29 #include <TList.h>
30 #include <TMathBase.h>
31 #include <TNamed.h>
32 #include <TObject.h>
33 #include <TPaveText.h>
34 #include <TROOT.h>
35 #include <TStopwatch.h>
36 #include <TString.h>
37 #include <TTask.h>
38 #include <utility>
39 
40 FairMonitor* FairMonitor::instance = nullptr;
41 
42 FairMonitor::FairMonitor()
43  : TNamed("FairMonitor", "Monitor for FairRoot")
44  , fRunMonitor(kFALSE)
45  , fDrawCanvas(kFALSE)
46  , fRunTime(0.)
47  , fRunMem(0.)
48  , fTimerMap()
49  , fMemoryMap()
50  , fHistList(new TList())
51  , fCanvas()
52  , fNoTaskRequired(0)
53  , fNoTaskCreated(0)
54  , fCurrentTask(0)
55  , fTaskRequired()
56  , fTaskCreated()
57  , fTaskCreatedTemp()
58  , fObjectMap()
59  , fTaskMap()
60  , fObjectPos()
61  , fTaskPos()
62 {}
63 
64 FairMonitor::~FairMonitor() {}
65 
67 {
68  if (!instance) {
69  instance = new FairMonitor();
70  }
71  return instance;
72 }
73 
74 void FairMonitor::StartTimer(const TTask* tTask, const char* identStr)
75 {
76  if (!fRunMonitor)
77  return;
78 
79  TString tempString = Form("timer_%p_%s_%s", tTask, tTask->GetName(), identStr);
80 
81  auto itt = fTimerMap.find(tempString);
82  if (itt == fTimerMap.end()) {
83  fTimerMap.insert(std::pair<TString, TStopwatch>(tempString, TStopwatch()));
84  itt = fTimerMap.find(tempString);
85  }
86  fTimerMap[tempString].Start();
87 }
88 
89 void FairMonitor::StopTimer(const TTask* tTask, const char* identStr)
90 {
91  if (!fRunMonitor)
92  return;
93 
94  TString tempString = Form("timer_%p_%s_%s", tTask, tTask->GetName(), identStr);
95 
96  auto itt = fTimerMap.find(tempString);
97  if (itt == fTimerMap.end()) {
98  LOG(info) << "FairMonitor::StopTimer() called without matching StartTimer()";
99  return;
100  }
101  // itt->second
102  fTimerMap[tempString].Stop();
103 
104  Double_t time = fTimerMap[tempString].RealTime();
105  RecordInfo(tTask, Form("%s_TIM", identStr), time);
106 
107  if (tempString.EndsWith("_EXEC"))
108  fRunTime += time;
109 }
110 
111 void FairMonitor::StartMemoryMonitor(const TTask* tTask, const char* identStr)
112 {
113  if (!fRunMonitor)
114  return;
115  FairSystemInfo sysInfo;
116  Int_t memoryAtStart = static_cast<Int_t>(sysInfo.GetCurrentMemory());
117 
118  TString memMon = Form("mem_%p_%s_%s", tTask, tTask->GetName(), identStr);
119 
120  auto iti = fMemoryMap.find(memMon);
121  if (iti == fMemoryMap.end())
122  fMemoryMap.insert(std::pair<TString, Int_t>(memMon, memoryAtStart));
123  else
124  iti->second = memoryAtStart;
125 }
126 
127 void FairMonitor::StopMemoryMonitor(const TTask* tTask, const char* identStr)
128 {
129  if (!fRunMonitor)
130  return;
131  FairSystemInfo sysInfo;
132  Int_t memoryAtEnd = static_cast<Int_t>(sysInfo.GetCurrentMemory());
133 
134  TString memMon = Form("mem_%p_%s_%s", tTask, tTask->GetName(), identStr);
135 
136  auto iti = fMemoryMap.find(memMon);
137  if (iti == fMemoryMap.end())
138  LOG(warn) << "FairMonitor::StopMemoryMonitor() Could not find corresponding entry for \"" << memMon.Data()
139  << "\".";
140  else {
141  RecordInfo(tTask, Form("%s_MEM", identStr), memoryAtEnd - iti->second);
142 
143  if (memMon.EndsWith("_EXEC"))
144  fRunMem += memoryAtEnd - iti->second;
145  }
146 }
147 
148 void FairMonitor::RecordInfo(const TTask* tTask, const char* identStr, Double_t value)
149 {
150  if (!fRunMonitor)
151  return;
152 
153  TString tempString = Form("hist_%p_%s_%s", tTask, tTask->GetName(), identStr);
154 
155  Int_t nofHists = fHistList->GetEntries();
156  Int_t ihist = 0;
157  for (ihist = 0; ihist < nofHists; ++ihist) {
158  if (!tempString.CompareTo(fHistList->At(ihist)->GetName())) {
159  break;
160  }
161  }
162  if (ihist == nofHists) {
163  TString titleString = Form("Histogram %s for %s", identStr, tTask->GetName());
164  fHistList->Add(new TH1F(tempString, titleString, 1000, 0, 1000));
165  }
166  TH1F* tempHist = (static_cast<TH1F*>((fHistList->At(ihist))));
167  Int_t nofEntries = tempHist->GetEntries();
168  if (nofEntries > tempHist->GetNbinsX())
169  tempHist->SetBins(tempHist->GetNbinsX() * 10, 0, tempHist->GetXaxis()->GetXmax() * 10);
170 
171  tempHist->SetBinContent(nofEntries + 1, value);
172 }
173 
174 void FairMonitor::RecordRegister(const char* name, const char* folderName, Bool_t toFile)
175 {
176  if (!fRunMonitor)
177  return;
178 
179  LOG(debug) << "*** FM::RecordRegister(" << name << ", " << folderName << (toFile ? ", kTRUE)" : ", kFALSE")
180  << " for task >>" << fCurrentTask << "<< (" << (fCurrentTask ? fCurrentTask->GetName() : "") << ")";
181  if (fCurrentTask == 0) {
182  fNoTaskCreated++;
183  LOG(info) << "*** FM::RecordRegister(" << name << ", " << folderName << (toFile ? ", kTRUE)" : ", kFALSE")
184  << " WITHOUT TASK";
185  return;
186  }
187  TString tempString = Form("%p_%s", fCurrentTask, fCurrentTask->GetName());
188  if (toFile) {
189  fTaskCreated.insert(std::pair<TString, TString>(tempString, TString(name)));
190  } else {
191  fTaskCreatedTemp.insert(std::pair<TString, TString>(tempString, TString(name)));
192  }
193 }
194 
195 void FairMonitor::RecordGetting(const char* name)
196 {
197  if (!fRunMonitor)
198  return;
199 
200  LOG(debug) << "*** FM::RecordGetting(" << name << ") for task >>" << fCurrentTask << "<< ("
201  << (fCurrentTask ? fCurrentTask->GetName() : "") << ")";
202 
203  if (fCurrentTask == 0) {
204  fNoTaskRequired++;
205  return;
206  }
207  TString tempString = Form("%p_%s", fCurrentTask, fCurrentTask->GetName());
208  fTaskRequired.insert(std::pair<TString, TString>(tempString, TString(name)));
209 }
210 
211 void FairMonitor::PrintTask(TTask* tempTask, Int_t taskLevel) const
212 {
213  if (!fRunMonitor)
214  return;
215 
216  Int_t nofHists = fHistList->GetEntries();
217  TString tempString = "";
218 
219  Double_t timInt = -1.;
220  Double_t timEnt = 0.;
221  Int_t memInt = -1.;
222 
223  char byteChar[4] = {'B', 'k', 'M', 'G'};
224 
225  for (Int_t ihist = 0; ihist < nofHists; ihist++) {
226  tempString = Form("%s", fHistList->At(ihist)->GetName());
227  if (tempString.Contains(Form("%p", tempTask))) {
228  if (tempString.Contains("EXEC_TIM")) {
229  timInt = (static_cast<TH1F*>(fHistList->At(ihist))->Integral());
230  timEnt = (static_cast<TH1F*>(fHistList->At(ihist))->GetEntries());
231  }
232  if (tempString.Contains("EXEC_MEM")) {
233  memInt = (static_cast<TH1F*>(fHistList->At(ihist))->Integral());
234  }
235  }
236  }
237  if (timInt < 0) {
238  LOG(warn) << "FairMonitor::PrintTask(), task \"" << tempTask->GetName() << "\" not found!";
239  } else {
240  // LOG(info) << "\"" << tempTask->GetName() << "\" --> TIME integral = " << timInt << " -- MEMORY integral =
241  // " << memInt;
242  TString printString = Form("%f", timInt / timEnt);
243  for (Int_t itemp = printString.Length(); itemp < 10; itemp++)
244  printString.Insert(0, ' ');
245  if (printString.Length() > 10)
246  printString.Remove(11, 100);
247  printString += " s/ev |";
248  Double_t timePerc = 100. * timInt / fRunTime;
249  if (timePerc < 100.)
250  printString += ' ';
251  if (timePerc < 10.)
252  printString += ' ';
253  tempString = Form("%f", timePerc);
254  for (Int_t itemp = tempString.Length(); itemp < 8; itemp++)
255  printString += ' ';
256  printString += tempString;
257  printString.Remove(27, 100);
258  Int_t byteIdent = 0;
259  while (memInt > 1024) {
260  memInt = memInt / 1024;
261  byteIdent++;
262  }
263  tempString = Form("%4d", memInt);
264  printString += " % ] ";
265  if (memInt < 0)
266  printString += "- - -";
267  else {
268  printString += tempString;
269  printString += byteChar[byteIdent];
270  }
271  printString += " \"";
272  for (Int_t ilev = 0; ilev < taskLevel; ilev++)
273  printString += " ";
274  printString += tempTask->GetName();
275  if (printString.Length() > 80)
276  printString.Remove(80, 100);
277  printString += "\"";
278  Int_t timeFrac = static_cast<Int_t>((timePerc / 100. * 30.));
279 
280  printString.Insert(37, "\033[0m");
281  switch (byteIdent) {
282  case 0:
283  printString.Insert(32, "\033[42m");
284  break;
285  case 1:
286  printString.Insert(32, "\033[43m");
287  break;
288  default:
289  printString.Insert(32, "\033[41m");
290  break;
291  }
292  printString.Insert(timeFrac, "\033[0m");
293  if (timePerc < 30)
294  printString.Insert(0, "[\033[42m");
295  else if (timePerc < 90)
296  printString.Insert(0, "[\033[43m");
297  else
298  printString.Insert(0, "[\033[41m");
299  LOG(INFO) << printString.Data();
300  }
301 
302  TList* subTaskList = tempTask->GetListOfTasks();
303  if (!subTaskList)
304  return;
305  for (Int_t itask = 0; itask < subTaskList->GetEntries(); itask++) {
306  TTask* subTask = static_cast<TTask*>(subTaskList->At(itask));
307  if (subTask)
308  PrintTask(subTask, taskLevel + 1);
309  }
310 }
311 
312 void FairMonitor::Print(Option_t*) const
313 {
314  if (!fRunMonitor) {
315  LOG(warn) << "FairMonitor was disabled. Nothing to print!";
316  return;
317  }
318 
319  LOG(info) << "- Total Run Time: " << fRunTime << " s ---------------------------------------------------------";
320  TTask* mainFairTask = static_cast<TTask*>((gROOT->GetListOfBrowsables()->FindObject("FairTaskList")));
321  if (mainFairTask)
322  PrintTask(mainFairTask, 0);
323  LOG(info) << "-------------------------------------------------------------------------------------";
324 }
325 
326 void FairMonitor::PrintTask(TString specString) const
327 {
328  if (!fRunMonitor) {
329  LOG(warn) << "FairMonitor was disabled. Nothing to print!";
330  return;
331  }
332 
333  TString unitString = "";
334 
335  Int_t nofHists = fHistList->GetEntries();
336  for (Int_t ihist = 0; ihist < nofHists; ihist++) {
337  TString histString = Form("%s", fHistList->At(ihist)->GetName());
338  if (histString.Contains(specString)) {
339  histString.Replace(0, histString.First('_') + 1, "");
340  histString.Replace(0, histString.First('_') + 1, "");
341  Double_t integral = (static_cast<TH1F*>((fHistList->At(ihist)))->Integral());
342  Double_t entries = (static_cast<TH1F*>((fHistList->At(ihist)))->GetEntries());
343  // Double_t valPent = integral/entries;
344  if (histString.EndsWith("_TIM"))
345  unitString = " s";
346  if (histString.EndsWith("_MEM"))
347  unitString = " B";
348  LOG(info) << histString.Data() << " >>>>> " << integral << unitString.Data() << " / " << entries
349  << " ent = " << integral / entries << unitString.Data() << " / ent";
350  }
351  }
352 }
353 
354 void FairMonitor::Draw(Option_t*)
355 {
356  if (!fRunMonitor) {
357  LOG(warn) << "FairMonitor was disabled. Nothing to print!";
358  return;
359  }
360 
361  typedef std::map<TString, Int_t>::iterator tiMapIter;
362  tiMapIter iti;
363 
364  TTask* mainFairTask = static_cast<TTask*>((gROOT->GetListOfBrowsables()->FindObject("FairTaskList")));
365  if (!mainFairTask)
366  return;
367 
368  GetTaskMap(mainFairTask);
369 
370  for (auto itb = fTaskRequired.begin(); itb != fTaskRequired.end(); ++itb) {
371  iti = fTaskMap.find(itb->first);
372  if (iti == fTaskMap.end())
373  fTaskMap.insert(std::pair<TString, Int_t>(itb->first, 0));
374  iti = fObjectMap.find(itb->second);
375  if (iti == fObjectMap.end())
376  fObjectMap.insert(std::pair<TString, Int_t>(itb->second, 0));
377  }
378  for (auto itb = fTaskCreated.begin(); itb != fTaskCreated.end(); ++itb) {
379  iti = fTaskMap.find(itb->first);
380  if (iti == fTaskMap.end())
381  fTaskMap.insert(std::pair<TString, Int_t>(itb->first, 0));
382  iti = fObjectMap.find(itb->second);
383  if (iti == fObjectMap.end())
384  fObjectMap.insert(std::pair<TString, Int_t>(itb->second, 0));
385  }
386  for (auto itb = fTaskCreatedTemp.begin(); itb != fTaskCreatedTemp.end(); ++itb) {
387  iti = fTaskMap.find(itb->first);
388  if (iti == fTaskMap.end())
389  fTaskMap.insert(std::pair<TString, Int_t>(itb->first, 0));
390  iti = fObjectMap.find(itb->second);
391  if (iti == fObjectMap.end())
392  fObjectMap.insert(std::pair<TString, Int_t>(itb->second, 0));
393  }
394 
395  AnalyzeObjectMap(mainFairTask);
396 
397  Int_t maxHierarchyNumber = 0;
398  for (iti = fObjectMap.begin(); iti != fObjectMap.end(); ++iti) {
399  if (maxHierarchyNumber < iti->second)
400  maxHierarchyNumber = iti->second;
401  }
402 
403  LOG(debug) << "Max hierarchy number is " << maxHierarchyNumber;
404 
405  fCanvas = new TCanvas("MonitorCanvas", "Fair Monitor", 10, 10, 960, 600);
406  fCanvas->cd();
407 
408  fCanvas->SetFillStyle(4000);
409  fCanvas->Range(0, 0, 960, 800);
410  fCanvas->SetFillColor(0);
411  fCanvas->SetBorderSize(0);
412  fCanvas->SetBorderMode(0);
413  fCanvas->SetFrameFillColor(0);
414 
415  for (Int_t ihier = 0; ihier < maxHierarchyNumber + 1; ++ihier) {
416  Int_t nofHier = 0;
417  for (iti = fTaskMap.begin(); iti != fTaskMap.end(); ++iti) {
418  if (iti->second == ihier) {
419  nofHier++;
420  }
421  }
422  LOG(debug) << "There are " << nofHier << " tasks on level " << ihier << ".";
423 
424  if (ihier == 0) {
425  Double_t iObj = 0.;
426  for (iti = fTaskMap.begin(); iti != fTaskMap.end(); ++iti) {
427  if (iti->second == ihier) {
428  std::pair<Double_t, Double_t> tempPos(50, 600 - iObj * 40);
429  fTaskPos.insert(std::pair<TString, std::pair<Double_t, Double_t>>(iti->first, tempPos));
430  iObj += 1.;
431  }
432  }
433  } else {
434  Int_t iObj = 0;
435  Int_t secLine = 0;
436  Int_t secLineEven = 0;
437  if (nofHier > 9) {
438  secLine = 1;
439  if (nofHier % 2 == 0)
440  secLineEven = 1;
441  }
442  Int_t startingPosition = 575 - nofHier / (1 + secLine) * 45 - secLine * (1 - secLineEven) * 45;
443 
444  Double_t topEdge =
445  800.
446  - 800. * (2. * static_cast<Double_t>(ihier) - 0.5) / (static_cast<Double_t>(2 * maxHierarchyNumber + 1))
447  + secLine * 15;
448  LOG(debug) << "for level " << ihier << " will put top edge at " << topEdge << ". "
449  << (secLineEven ? "Two lines" : "One line") << (secLineEven ? " with offset" : "");
450  for (iti = fTaskMap.begin(); iti != fTaskMap.end(); ++iti) {
451  if (iti->second == ihier) {
452  std::pair<Double_t, Double_t> tempPos(startingPosition + iObj * 90 / (1 + secLine)
453  - secLineEven * (iObj % 2) * 45,
454  topEdge - 15 - secLine * (iObj % 2) * 35);
455  fTaskPos.insert(std::pair<TString, std::pair<Double_t, Double_t>>(iti->first, tempPos));
456  iObj += 1;
457  }
458  }
459  }
460 
461  nofHier = 0;
462  for (iti = fObjectMap.begin(); iti != fObjectMap.end(); ++iti) {
463  if (TMath::Abs(iti->second) == ihier) {
464  nofHier++;
465  }
466  }
467  LOG(debug) << "There are " << nofHier << " objects on level " << ihier << ".";
468 
469  Int_t iObj = 0;
470  Int_t secLine = 0;
471  Int_t secLineEven = 0;
472  if (nofHier > 9) {
473  secLine = 1;
474  if (nofHier % 2 == 0)
475  secLineEven = 1;
476  }
477  Int_t startingPosition = 575 - nofHier / (1 + secLine) * 45 - secLine * (1 - secLineEven) * 45;
478 
479  Double_t topEdge =
480  800.
481  - 800. * (2. * static_cast<Double_t>(ihier) + 0.5) / (static_cast<Double_t>(2 * maxHierarchyNumber + 1))
482  + secLine * 15;
483  LOG(debug) << "for level " << ihier << " will put top edge at " << topEdge << ". "
484  << (secLineEven ? "Two lines" : "One line") << (secLineEven ? " with offset" : "");
485  for (iti = fObjectMap.begin(); iti != fObjectMap.end(); ++iti) {
486  if (TMath::Abs(iti->second) == ihier) {
487  std::pair<Double_t, Double_t> tempPos(startingPosition + iObj * 90 / (1 + secLine)
488  - secLineEven * (iObj % 2) * 45,
489  topEdge - 15 - secLine * (iObj % 2) * 35);
490  fObjectPos.insert(std::pair<TString, std::pair<Double_t, Double_t>>(iti->first, tempPos));
491  iObj += 1;
492  }
493  }
494  }
495 
496  typedef std::map<TString, std::pair<Double_t, Double_t>>::iterator tddMapIter;
497  tddMapIter itt;
498  tddMapIter ito;
499 
500  for (auto itb = fTaskRequired.begin(); itb != fTaskRequired.end(); ++itb) {
501  itt = fTaskPos.find(itb->first);
502  ito = fObjectPos.find(itb->second);
503  if (itt != fTaskPos.end() && ito != fObjectPos.end()) {
504  std::pair<Double_t, Double_t> taskPos = itt->second;
505  std::pair<Double_t, Double_t> objectPos = ito->second;
506  TArrow* tempArrow =
507  new TArrow(objectPos.first, objectPos.second - 15, taskPos.first, taskPos.second + 15, 0.01, "|>");
508  tempArrow->Draw();
509  }
510  }
511  for (auto itb = fTaskCreated.begin(); itb != fTaskCreated.end(); ++itb) {
512  itt = fTaskPos.find(itb->first);
513  ito = fObjectPos.find(itb->second);
514  if (itt != fTaskPos.end() && ito != fObjectPos.end()) {
515  std::pair<Double_t, Double_t> taskPos = itt->second;
516  std::pair<Double_t, Double_t> objectPos = ito->second;
517  TArrow* tempArrow =
518  new TArrow(taskPos.first, taskPos.second - 15, objectPos.first, objectPos.second + 15, 0.01, "|>");
519  tempArrow->Draw();
520  }
521  }
522  for (auto itb = fTaskCreatedTemp.begin(); itb != fTaskCreatedTemp.end(); ++itb) {
523  itt = fTaskPos.find(itb->first);
524  ito = fObjectPos.find(itb->second);
525  if (itt != fTaskPos.end() && ito != fObjectPos.end()) {
526  std::pair<Double_t, Double_t> taskPos = itt->second;
527  std::pair<Double_t, Double_t> objectPos = ito->second;
528  TArrow* tempArrow =
529  new TArrow(taskPos.first, taskPos.second - 15, objectPos.first, objectPos.second + 15, 0.01, "|>");
530  tempArrow->Draw();
531  }
532  }
533 
534  for (itt = fTaskPos.begin(); itt != fTaskPos.end(); ++itt) {
535  std::pair<Double_t, Double_t> taskPos = itt->second;
536  TBox* bkgBox = new TBox(taskPos.first - 40, taskPos.second - 15, taskPos.first + 40, taskPos.second + 15);
537  bkgBox->SetFillColor(kGreen - 9);
538  bkgBox->Draw();
539 
540  TString tempString = Form("hist_%s_%s", itt->first.Data(), "EXEC_TIM");
541  Int_t nofHists = fHistList->GetEntries();
542  for (Int_t ihist = 0; ihist < nofHists; ++ihist) {
543  if (!tempString.CompareTo(fHistList->At(ihist)->GetName())) {
544  Double_t timeInt = (static_cast<TH1F*>(fHistList->At(ihist))->Integral());
545  Double_t timeFrac = 80. * timeInt / fRunTime;
546  TBox* barBox = new TBox(
547  taskPos.first - 40, taskPos.second - 15, taskPos.first - 40 + timeFrac, taskPos.second + 15);
548  barBox->SetFillColor(kRed);
549  barBox->Draw();
550  }
551  }
552 
553  tempString = itt->first;
554  tempString.Replace(0, tempString.First('_') + 1, "");
555  if (tempString.Length() > 16)
556  tempString.Replace(16, tempString.Length(), "");
557  TLatex* taskText = new TLatex(taskPos.first, taskPos.second, tempString.Data());
558  taskText->SetTextAlign(22);
559  taskText->SetTextSize(0.015);
560  taskText->Draw();
561  }
562 
563  for (ito = fObjectPos.begin(); ito != fObjectPos.end(); ++ito) {
564  std::pair<Double_t, Double_t> objectPos = ito->second;
565  iti = fObjectMap.find(ito->first);
566  TPaveText* paveText = new TPaveText(
567  objectPos.first - 40, objectPos.second - 15, objectPos.first + 40, objectPos.second + 15, "nb");
568  paveText->SetFillColor(kMagenta - 9);
569  if (iti != fObjectMap.end())
570  if (iti->second < 0)
571  paveText->SetFillColor(kGray);
572  paveText->SetShadowColor(0);
573  paveText->SetTextSize(0.015);
574  TString tempString = ito->first;
575  if (tempString.Length() > 16)
576  tempString.Replace(16, tempString.Length(), "");
577  paveText->AddText(tempString);
578  paveText->Draw();
579  }
580 }
581 
582 void FairMonitor::DrawHist(TString specString)
583 {
584  if (!fRunMonitor) {
585  LOG(warn) << "FairMonitor was disabled. Nothing to draw!";
586  return;
587  }
588 
589  Int_t nofHists = fHistList->GetEntries();
590  Int_t nofDraws = 0;
591  TString drawStr = "P";
592 
593  TLegend* tempLeg = new TLegend(0.0, 0.5, 0.5, 0.9);
594 
595  Double_t histMax = 0.;
596  for (Int_t ihist = 0; ihist < nofHists; ihist++) {
597  TString histString = Form("%s", fHistList->At(ihist)->GetName());
598  if (histString.Contains(specString)) {
599  TH1F* tempHist = (static_cast<TH1F*>((fHistList->At(ihist))));
600  if (histMax < tempHist->GetMaximum())
601  histMax = tempHist->GetMaximum();
602  }
603  }
604 
605  for (Int_t ihist = 0; ihist < nofHists; ihist++) {
606  TString histString = Form("%s", fHistList->At(ihist)->GetName());
607  if (histString.Contains(specString)) {
608  TH1F* tempHist = (static_cast<TH1F*>((fHistList->At(ihist))));
609  if (tempHist->GetXaxis()->GetXmax() > tempHist->GetEntries())
610  tempHist->SetBins(tempHist->GetEntries(), 0, tempHist->GetEntries());
611  tempHist->SetMarkerColor(nofDraws % 6 + 2);
612  tempHist->SetMarkerStyle(nofDraws % 4 + 20);
613  tempHist->SetAxisRange(0., histMax * 1.1, "Y");
614  tempHist->Draw(drawStr);
615  TString tempName = tempHist->GetName();
616  tempName.Remove(0, tempName.First('_') + 1);
617  tempName.Remove(0, tempName.First('_') + 1);
618  tempLeg->AddEntry(tempHist, tempName, "p");
619  nofDraws++;
620  drawStr = "Psame";
621  }
622  }
623  if (nofDraws > 1)
624  tempLeg->Draw();
625 }
626 
627 void FairMonitor::StoreHistograms(TFile* sinkFile)
628 {
629  if (!fRunMonitor) {
630  return;
631  }
632  Bool_t wasBatch = gROOT->IsBatch();
633  if (!fDrawCanvas && !wasBatch) // default: switching to batch mode if draw canvas disabled
634  gROOT->SetBatch(kTRUE);
635  this->Draw();
636  if (!fDrawCanvas && !wasBatch) // default: switching to batch mode if draw canvas disabled
637  gROOT->SetBatch(kFALSE);
638 
639  TFile* histoFile = sinkFile;
640  if (fOutputFileName.Length() > 1 && fOutputFileName != sinkFile->GetName()) {
641  histoFile = TFile::Open(fOutputFileName, "recreate");
642  }
643 
644  histoFile->mkdir("MonitorResults");
645  histoFile->cd("MonitorResults");
646  TIter next(fHistList);
647  while (TH1* thist = (static_cast<TH1*>(next()))) {
648  thist->SetBins(thist->GetEntries(), 0, thist->GetEntries());
649  thist->Write();
650  }
651  if (fCanvas) {
652  fCanvas->Write();
653  }
654 
655  if (histoFile != sinkFile) {
656  histoFile->Close();
657  delete histoFile;
658  }
659 }
660 
661 // Private function to fill the map of the tasks
662 void FairMonitor::GetTaskMap(TTask* tempTask)
663 {
664  TString tempString = Form("%p_%s", tempTask, tempTask->GetName());
665  fTaskMap.insert(std::pair<TString, Int_t>(tempString, 0));
666 
667  TList* subTaskList = tempTask->GetListOfTasks();
668  if (!subTaskList)
669  return;
670  for (Int_t itask = 0; itask < subTaskList->GetEntries(); itask++) {
671  TTask* subTask = static_cast<TTask*>(subTaskList->At(itask));
672  if (subTask)
673  GetTaskMap(subTask);
674  }
675 }
676 
677 // Private function to analyze the object list:
678 // assign to each of the tasks and objects (TClonesArray, pressumably)
679 // one number that places them in a hierarchy
680 // objects number 0 are the ones in the input array
681 // tasks get number from the abs(highest number object they read) increased by one
682 // tasks number i create objects number i (if they go to output file) or
683 // tasks number i create objects number -i (if they are not persistent)
684 void FairMonitor::AnalyzeObjectMap(TTask* tempTask)
685 {
686  TString tempString = Form("%p_%s", tempTask, tempTask->GetName());
687 
688  Int_t hierarchyNumber = 0;
689 
690  typedef std::map<TString, Int_t>::iterator tiMapIter;
691  tiMapIter iti;
692 
693  LOG(debug) << "TASK \"" << tempTask->GetName() << "\" NEEDS:";
694  for (auto itb = fTaskRequired.begin(); itb != fTaskRequired.end(); ++itb) {
695  if (itb->first != tempString)
696  continue;
697  LOG(debug) << " \"" << itb->second.Data() << "\"";
698  iti = fObjectMap.find(itb->second);
699  if (iti == fObjectMap.end())
700  continue;
701  if (hierarchyNumber <= TMath::Abs(iti->second))
702  hierarchyNumber = TMath::Abs(iti->second) + 1;
703  }
704 
705  // hierarchyNumber++;
706 
707  LOG(debug) << "WILL GET hierarchyNumber = " << hierarchyNumber;
708 
709  iti = fTaskMap.find(tempString);
710  if (iti != fTaskMap.end())
711  iti->second = hierarchyNumber;
712 
713  LOG(debug) << " \"" << tempTask->GetName() << "\" CREATES:";
714  for (auto itb = fTaskCreated.begin(); itb != fTaskCreated.end(); ++itb) {
715  if (itb->first != tempString)
716  continue;
717  LOG(debug) << " + \"" << itb->second.Data() << "\"";
718  iti = fObjectMap.find(itb->second);
719  if (iti == fObjectMap.end())
720  continue;
721  iti->second = hierarchyNumber;
722  }
723 
724  for (auto itb = fTaskCreatedTemp.begin(); itb != fTaskCreatedTemp.end(); ++itb) {
725  if (itb->first != tempString)
726  continue;
727  LOG(debug) << " - \"" << itb->second.Data() << "\"";
728  iti = fObjectMap.find(itb->second);
729  if (iti == fObjectMap.end())
730  continue;
731  iti->second = -hierarchyNumber;
732  }
733 
734  TList* subTaskList = tempTask->GetListOfTasks();
735  if (!subTaskList)
736  return;
737  for (Int_t itask = 0; itask < subTaskList->GetEntries(); itask++) {
738  TTask* subTask = static_cast<TTask*>(subTaskList->At(itask));
739  if (subTask)
740  AnalyzeObjectMap(subTask);
741  }
742 }
743 
void RecordInfo(const TTask *tTask, const char *identStr, Double_t value)
virtual void Print(Option_t *option="") const
void DrawHist(TString specString)
void StartMemoryMonitor(const TTask *tTask, const char *identStr)
void StopTimer(const TTask *tTask, const char *identStr)
Definition: FairMonitor.cxx:89
static FairMonitor * GetMonitor()
Definition: FairMonitor.cxx:66
ClassImp(FairEventBuilder)
void RecordGetting(const char *name)
void StopMemoryMonitor(const TTask *tTask, const char *identStr)
void RecordRegister(const char *name, const char *folderName, Bool_t toFile)
void StartTimer(const TTask *tTask, const char *identStr)
Definition: FairMonitor.cxx:74
virtual void Draw(Option_t *option="")
size_t GetCurrentMemory()
void StoreHistograms(TFile *sinkFile)
void PrintTask(TString specString) const