Jpp  18.6.0-rc.1
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JGizmoToolkit.hh
Go to the documentation of this file.
1 #ifndef __JGIZMOTOOLKIT__
2 #define __JGIZMOTOOLKIT__
3 
4 #include <string>
5 #include <sstream>
6 #include <map>
7 #include <cmath>
8 #include <memory>
9 
10 #pragma GCC diagnostic push
11 #pragma GCC diagnostic ignored "-Wall"
12 #include "TError.h"
13 #include "TFile.h"
14 #include "TClass.h"
15 #include "TObject.h"
16 #include "TKey.h"
17 #include "TH1.h"
18 #include "TH2.h"
19 #include "TH3.h"
20 #include "TGraph.h"
21 #include "TGraphErrors.h"
22 #include "TGraph2D.h"
23 #include "TGraph2DErrors.h"
24 #include "TMultiGraph.h"
25 #include "TLine.h"
26 #include "TEllipse.h"
27 #include "TString.h"
28 #include "TRegexp.h"
29 #include "TFormula.h"
30 #include "TF1.h"
31 #include "TF2.h"
32 #include "TList.h"
33 #include "TIterator.h"
34 #include "TMethod.h"
35 #include "TMethodCall.h"
36 #include "TAxis.h"
37 #include "TMath.h"
38 #pragma GCC diagnostic pop
39 
40 #include "JLang/JException.hh"
41 #include "JGizmo/JRootObjectID.hh"
42 
44 
45 
46 /**
47  * \author mdejong
48  */
49 
50 namespace JGIZMO {}
51 namespace JPP { using namespace JGIZMO; }
52 
53 /**
54  * Auxiliary applications for use of ROOT and more.
55  */
56 namespace JGIZMO {
57 
58  using JLANG::JParseError;
62 
63 
64  /**
65  * Auxiliary data structure for JOpera1D.cc and JOpera2D.cc applications.
66  */
67  struct JOpera {
68  //
69  // Histogram name.
70  //
71  static const char* const SAME_AS_OPERATION() { return "%"; } //!< Set name of output histogram to name of operation
72  static const char* const SAME_AS_INPUT() { return "="; } //!< Set name of output histogram to name of input histogram
73 
74  //
75  // Histogram operations.
76  //
77  static const char* const Add() { return "Add"; } //!< ROOT TH1::Add
78  static const char* const add() { return "add"; } //!< Add contents with lookup bin in second histogram
79  static const char* const Subtract() { return "Subtract"; } //!< ROOT TH1::Subtract
80  static const char* const subtract() { return "subtract"; } //!< Subtract contents with lookup bin in second histogram
81  static const char* const Multiply() { return "Multiply"; } //!< ROOT TH1::Multiply
82  static const char* const multiply() { return "multiply"; } //!< Multiply contents with lookup bin in second histogram
83  static const char* const Divide() { return "Divide"; } //!< ROOT TH1::Divide
84  static const char* const divide() { return "divide"; } //!< Divide contents with lookup bin in second histogram
85  static const char* const efficiency() { return "efficiency"; } //!< Divide contents and multiply errors with inefficiency
86  static const char* const stdev() { return "stdev"; } //!< Set contents to standard deviation
87  static const char* const sqrt() { return "sqrt"; } //!< Set contents to signed difference between squares
88  static const char* const Replace() { return "Replace"; } //!< Set contents to associated function
89  };
90 
91 
92  /**
93  * Time stamp of earliest UTC time.
94  */
95  static const char* const TIMESTAMP = "#splitline{}{#splitline{%d:%m:%y}{ %H:%M}}%F1970-01-01 00:00:00";
96 
97 
98  /**
99  * Get TFile pointer corresponding to give file name.
100  *
101  * The TFile pointer of an already opened file is recovered, else a new file is opened.\n
102  * Note that the closure of the opened files should be done by the caller of this method.
103  *
104  * \param file_name file name
105  * \param option TFile::Open option
106  * \return pointer to TFile
107  */
108  inline TFile* getFile(const std::string& file_name, const std::string& option = "exist")
109  {
110  using namespace std;
111 
112  gErrorIgnoreLevel = kError;
113 
114  static map<string, TFile*> zmap;
115 
116  map<string, TFile*>::iterator i = zmap.find(file_name);
117 
118  if (i == zmap.end() || i->second == NULL || !i->second->IsOpen()) {
119 
120  TFile* file = TFile::Open(file_name.c_str(), option.c_str());
121 
122  zmap[file_name] = file;
123 
124  return file;
125 
126  } else {
127 
128  return i->second;
129  }
130  }
131 
132 
133  /**
134  * Get TDirectory pointer.
135  *
136  * The TFile pointer of an already opened file is recovered, else a new file is opened.
137  *
138  * \param id identifier
139  * \return pointer to TDirectory
140  */
141  inline TDirectory* getDirectory(const JRootObjectID& id)
142  {
143  TFile* in = getFile(id.getFilename().c_str(), "exist");
144 
145  if (in == NULL || !in->IsOpen()) {
146  return NULL;
147  }
148 
149  if (id.getDirectory() != "")
150  return in->GetDirectory(id.getDirectory());
151  else
152  return in;
153  }
154 
155 
156  /**
157  * Get first TObject with given identifier.
158  *
159  * \param id identifier
160  * \return pointer to TObject (or NULL)
161  */
162  inline TObject* getObject(const JRootObjectID& id)
163  {
164  TDirectory* dir = getDirectory(id);
165 
166  if (dir != NULL) {
167 
168  const TRegexp regexp(id.getObjectName());
169 
170  TIter iter(dir->GetListOfKeys());
171 
172  for (TKey* key; (key = (TKey*) iter.Next()) != NULL; ) {
173 
174  const TString tag(key->GetName());
175 
176  // option match
177 
178  if (tag.Index(regexp) != -1) {
179  return key->ReadObj();
180  }
181  }
182  }
183 
184  return NULL;
185  }
186 
187 
188  /**
189  * Get drawing option of TH1.
190  *
191  * \param object pointer to TObject
192  * \return true if TH1 looks like a line; else false
193  */
194  inline bool isTAttLine(const TObject* object)
195  {
196  {
197  const TH1* h1 = dynamic_cast<const TH1*>(object);
198 
199  if (h1 != NULL) {
200 
201  if (h1->GetSumw2N()) {
202  for (Int_t i = 1; i <= h1->GetNbinsX(); ++i) {
203  if (h1->GetBinError(i) != 0.0) {
204  return false;
205  }
206  }
207  }
208  }
209  }
210  {
211  const TGraphErrors* g1 = dynamic_cast<const TGraphErrors*>(object);
212 
213  if (g1 != NULL) {
214 
215  for (Int_t i = 0; i != g1->GetN(); ++i) {
216  if (g1->GetEY()[i] != 0.0) {
217  return false;
218  }
219  }
220 
221  return g1->GetN() > 1;
222  }
223  }
224  {
225  const TGraph* g1 = dynamic_cast<const TGraph*>(object);
226 
227  if (g1 != NULL) {
228  return g1->GetN() > 1;
229  }
230  }
231 
232  return true;
233  }
234 
235 
236  /**
237  * Get result of given textual formula.
238  *
239  * The formula may contain names of member methods of the object pointed to.\n
240  * These methods should return a value that is compatible with <tt>Double_t</tt> and could have arguments.\n
241  * For example:
242  * <pre>
243  * getResult("1.0 / GetEntries", TH1*);
244  * </pre>
245  *
246  * \param text text
247  * \param object pointer to object
248  * \return value
249  */
250  inline Double_t getResult(const TString& text, TObject* object = NULL)
251  {
252  TString buffer(text);
253 
254  if (object != NULL) {
255 
256  TClass* p = TClass::GetClass(object->ClassName());
257 
258  if (p != NULL) {
259 
260  for ( ; ; ) {
261 
262  TMethod* method = NULL;
263 
264  for (std::unique_ptr<TIterator> iter(p->GetListOfAllPublicMethods()->MakeIterator()); TMethod* p = (TMethod*) iter->Next(); ) {
265  if (buffer.Index(p->GetName()) != -1) {
266  if (method == NULL || strlen(p->GetName()) > strlen(method->GetName())) {
267  method = p;
268  }
269  }
270  }
271 
272  if (method == NULL) {
273  break;
274  }
275 
276  for (Ssiz_t index; (index = buffer.Index(method->GetName())) != -1; ) {
277 
278  const TRegexp fp(" *([^)]*)"); // function call
279 
280  Ssiz_t len;
281  Ssiz_t pos = buffer.Index(fp, &len, index);
282 
283  Double_t value;
284 
285  if (pos == -1 || pos != index + (Ssiz_t) strlen(method->GetName())) {
286 
287  TMethodCall(p, method->GetName(), NULL).Execute(object, value);
288 
289  len = strlen(method->GetName());
290 
291  } else {
292 
293  TMethodCall(p, method->GetName(), NULL).Execute(object, TString(buffer(pos + 1, len - 2)), value);
294 
295  len += strlen(method->GetName());
296  }
297 
298  buffer.Replace(index, len, TString::Format("%20.10e", value));
299  }
300  }
301  }
302  }
303 
304  return TFormula("/tmp", buffer.Data()).Eval(0.0);
305  }
306 
307 
308  /**
309  * Get result of given textual formula.
310  *
311  * The formula may contain names of member methods of the object pointed to.\n
312  * These methods should return a value that is compatible with <tt>Double_t</tt> and could have arguments.\n
313  * For example:
314  * <pre>
315  * getResult("1.0 / GetEntries", TH1*);
316  * </pre>
317  *
318  * \param text text
319  * \param object pointer to object
320  * \return value
321  */
322  inline Double_t getResult(const std::string& text, TObject* object = NULL)
323  {
324  return getResult(TString(text.c_str()), object);
325  }
326 
327 
328  /**
329  * Get parameter number from text string.
330  *
331  * The number corresponds to the value <tt>[0-9]*</tt> in the expression <tt>"p[0-9]* = .."</tt>.
332  *
333  * \param text text
334  * \return parameter number
335  */
336  inline int getParameter(const std::string& text)
337  {
338  const char* regexp("p[0-9]* *=");
339 
340  TString buffer(text.c_str());
341 
342  buffer = buffer(TRegexp(regexp));
343  buffer = buffer(1, buffer.Length() - 2);
344 
345  if (!buffer.IsDigit()) {
346  THROW(JParseError, "Text is not a number " << text << ' ' << regexp);
347  }
348 
349  return buffer.Atoi();
350  }
351 
352 
353  /**
354  * Get parameter value from text string.
355  *
356  * The formula may contain names of member methods of the object pointed to.\n
357  * These methods should return a value that is compatible with <tt>Double_t</tt> and could have arguments.\n
358  * For example:
359  * <pre>
360  * getValue("p[..] = 2 * GetMaximum", TH1*);
361  * </pre>
362  *
363  * \param text text
364  * \param object pointer to object
365  * \return value
366  */
367  inline Double_t getValue(const std::string& text, TObject* object = NULL)
368  {
369  const char* regexp("=.*");
370 
371  TString buffer(text.c_str());
372 
373  buffer = buffer(TRegexp(regexp));
374  buffer = buffer(1, buffer.Length() - 1);
375 
376  return getResult(std::string(buffer), object);
377  }
378 
379 
380  /**
381  * Get parameter value from text string.
382  *
383  * The formula may contain names of member methods of the object pointed to.\n
384  * These methods should return a value that is compatible with <tt>Double_t</tt> and could have arguments.\n
385  * For example:
386  * <pre>
387  * getValue("p[..] = 1.0 2.0 3.0", 1);
388  * </pre>
389  * will return <tt>2.0</tt>.
390  *
391  * \param text text
392  * \param index index
393  * \return value
394  */
395  inline Double_t getValue(const std::string& text, const int index)
396  {
397  using namespace std;
398 
399  const char* regexp("=.*");
400 
401  TString buffer(text.c_str());
402 
403  buffer = buffer(TRegexp(regexp));
404  buffer = buffer(1, buffer.Length() - 1);
405 
406 
407  istringstream is((std::string(buffer)));
408 
409  Double_t value;
410 
411  for (int i = index; is >> value && i > 0; --i) {}
412 
413  if (is)
414  return value;
415  else
416  THROW(JParseError, "Text does not contain a number at given position " << buffer << ' ' << index);
417  }
418 
419 
420  /**
421  * Make histogram axis logarithmic (e.g.\ after using <tt>log10()</tt>).
422  *
423  * \param axis axis
424  */
425  inline void setLogarithmic(TAxis* axis)
426  {
427  if (axis != NULL) {
428 
429  const Int_t first = axis->GetFirst();
430  const Int_t last = axis->GetLast();
431 
432  const Double_t xmin = axis->GetBinLowEdge(first);
433  const Double_t xmax = axis->GetBinLowEdge(last) + axis->GetBinWidth(last);
434 
435  const Int_t N = axis->GetNbins();
436  Double_t buffer[N+1];
437 
438  buffer[0] = TMath::Power(10.0, axis->GetBinLowEdge(1));
439 
440  for (Int_t i = 1; i <= N; ++i) {
441  buffer[i] = TMath::Power(10.0, axis->GetBinLowEdge(i) + axis->GetBinWidth(i));
442  }
443 
444  axis->Set(N, buffer);
445 
446  if (axis->TestBit(TAxis::kAxisRange)) {
447  axis->SetRangeUser(TMath::Power(10.0, xmin), TMath::Power(10.0, xmax));
448  }
449  }
450  }
451 
452 
453  /**
454  * Make given parameter in formula logarithmic (e.g.\ after using <tt>log10()</tt>).
455  *
456  * \param formula formula
457  * \param parameter parameter
458  * \return formula
459  */
460  inline TString getLogarithmic(const TString& formula, const char parameter)
461  {
462  const TRegexp regexp[] = {
463  TString("^") + TString(parameter) + TString("[^a-zA-Z_]"), // parameter at start of line
464  TString("[^a-zA-Z_]") + TString(parameter) + TString("[^a-zA-Z_]"), // parameter in middle of line
465  TString("[^a-zA-Z_]") + TString(parameter) + TString("$") // parameter at end of line
466  };
467 
468  const TString replacement = TString("log10(") + TString(parameter) + TString(")");
469 
470  TString buffer(formula);
471 
472  if (buffer.Length() == 1 && buffer[0] == parameter) {
473 
474  buffer = replacement;
475 
476  } else {
477 
478  for (Ssiz_t pos = 0, i; pos < buffer.Length(); pos += replacement.Length()) {
479  if ((i = buffer.Index(regexp[0], pos)) != -1) { buffer.Replace((pos = i + 0), 1, replacement); }
480  else if ((i = buffer.Index(regexp[1], pos)) != -1) { buffer.Replace((pos = i + 1), 1, replacement); }
481  else if ((i = buffer.Index(regexp[2], pos)) != -1) { buffer.Replace((pos = i + 1), 1, replacement); }
482  else { break; }
483  }
484  }
485 
486  return buffer;
487  }
488 
489 
490  /**
491  * Copy function parameters.
492  *
493  * \param from function
494  * \param to function
495  */
496  inline void copy(const TF1& from, TF1& to)
497  {
498  static_cast<TAttLine&> (to) = static_cast<const TAttLine&> (from);
499  static_cast<TAttFill&> (to) = static_cast<const TAttFill&> (from);
500  static_cast<TAttMarker&>(to) = static_cast<const TAttMarker&>(from);
501 
502  to.SetParameters(from.GetParameters());
503 
504  to.SetNpx(from.GetNpx());
505  }
506 
507 
508  /**
509  * Copy function parameters.
510  *
511  * \param from function
512  * \param to function
513  */
514  inline void copy(const TF2& from, TF2& to)
515  {
516  copy(static_cast<const TF1&>(from), static_cast<TF1&>(to));
517 
518  to.SetNpy(from.GetNpy());
519  }
520 
521 
522  /**
523  * Make x-axis of objects in list logarithmic (e.g.\ after using <tt>log10()</tt>).
524  *
525  * \param list list
526  */
527  template<class T>
528  inline void setLogarithmicX(TList* list);
529 
530 
531  /**
532  * Make y-axis of objects in list logarithmic (e.g.\ after using <tt>log10()</tt>).
533  *
534  * \param list list
535  */
536  template<class T>
537  inline void setLogarithmicY(TList* list);
538 
539 
540  /**
541  * Make x-axis of given function logarithmic (e.g.\ after using <tt>log10()</tt>).
542  *
543  * \param f1 function
544  */
545  inline void setLogarithmicX(TF1* f1)
546  {
547  if (f1 != NULL) {
548 
549  TF1 fn(f1->GetName(), getLogarithmic(f1->GetExpFormula(), 'x'));
550 
551  copy(*f1, fn);
552 
553  fn.SetRange(f1->GetXmin(),
554  f1->GetXmax());
555 
556  *f1 = fn;
557 
558  f1->SetRange(TMath::Power(10.0,f1->GetXmin()),
559  TMath::Power(10.0,f1->GetXmax()));
560  }
561  }
562 
563 
564  /**
565  * Make y-axis of given function logarithmic (e.g.\ after using <tt>log10()</tt>).
566  *
567  * \param f1 function
568  */
569  inline void setLogarithmicY(TF1* f1)
570  {
571  if (f1 != NULL) {
572 
573  TString buffer = f1->GetExpFormula();
574 
575  buffer = "pow(10.0, " + buffer + ")";
576 
577  TF1 fn(f1->GetName(), buffer);
578 
579  copy(*f1, fn);
580 
581  *f1 = fn;
582  }
583  }
584 
585 
586  /**
587  * Make x-axis of given function logarithmic (e.g.\ after using <tt>log10()</tt>).
588  *
589  * \param f2 function
590  */
591  inline void setLogarithmicX(TF2* f2)
592  {
593  if (f2 != NULL) {
594 
595  TF2 fn(f2->GetName(), getLogarithmic(f2->GetExpFormula(), 'x'));
596 
597  copy(*f2, fn);
598 
599  fn.SetRange(f2->GetXmin(),
600  f2->GetYmin(),
601  f2->GetXmax(),
602  f2->GetYmax());
603 
604  *f2 = fn;
605 
606  f2->SetRange(TMath::Power(10.0,f2->GetXmin()),
607  f2->GetYmin(),
608  TMath::Power(10.0,f2->GetXmax()),
609  f2->GetYmax());
610  }
611  }
612 
613 
614  /**
615  * Make y-axis of given function logarithmic (e.g.\ after using <tt>log10()</tt>).
616  *
617  * \param f2 function
618  */
619  inline void setLogarithmicY(TF2* f2)
620  {
621  if (f2 != NULL) {
622 
623  TF2 fn(f2->GetName(), getLogarithmic(f2->GetExpFormula(), 'y'));
624 
625  copy(*f2, fn);
626 
627  fn.SetRange(f2->GetXmin(),
628  f2->GetYmin(),
629  f2->GetXmax(),
630  f2->GetYmax());
631 
632  *f2 = fn;
633 
634  f2->SetRange(f2->GetXmin(),
635  TMath::Power(10.0,f2->GetYmin()),
636  f2->GetXmax(),
637  TMath::Power(10.0,f2->GetYmax()));
638  }
639  }
640 
641 
642  /**
643  * Make x-axis and associated functions of given histogram logarithmic (e.g.\ after using <tt>log10()</tt>).
644  *
645  * \param h1 histogram
646  */
647  inline void setLogarithmicX(TH1* h1)
648  {
649  if (h1 != NULL) {
650 
651  setLogarithmicX<TF1>(h1->GetListOfFunctions());
652 
653  setLogarithmic(h1->GetXaxis());
654  }
655  }
656 
657 
658  /**
659  * Make x-axis and associated functions of given histogram logarithmic (e.g.\ after using <tt>log10()</tt>).
660  *
661  * \param h2 histogram
662  */
663  inline void setLogarithmicX(TH2* h2)
664  {
665  using namespace std;
666 
667  if (h2 != NULL) {
668 
669  setLogarithmicX<TF2>(h2->GetListOfFunctions());
670 
671  setLogarithmic(h2->GetXaxis());
672  }
673  }
674 
675 
676  /**
677  * Make y-axis and associated functions of given histogram logarithmic (e.g.\ after using <tt>log10()</tt>).
678  *
679  * \param h2 histogram
680  */
681  inline void setLogarithmicY(TH2* h2)
682  {
683  if (h2 != NULL) {
684 
685  setLogarithmicY<TF2>(h2->GetListOfFunctions());
686 
687  setLogarithmic(h2->GetYaxis());
688  }
689  }
690 
691 
692  /**
693  * Make x-axis and associated functions of given graph logarithmic (e.g.\ after using <tt>log10()</tt>).
694  *
695  * \param g1 graph
696  */
697  inline void setLogarithmicX(TGraph* g1)
698  {
699  if (g1 != NULL) {
700 
701  setLogarithmicX<TF1>(g1->GetListOfFunctions());
702 
703  for (Int_t i = 0; i != g1->GetN(); ++i) {
704  g1->GetX()[i] = pow(10.0, g1->GetX()[i]);
705  }
706  }
707  }
708 
709 
710  /**
711  * Make y-axis and associated functions of given graph logarithmic (e.g.\ after using <tt>log10()</tt>).
712  *
713  * \param g1 graph
714  */
715  inline void setLogarithmicY(TGraph* g1)
716  {
717  if (g1 != NULL) {
718 
719  setLogarithmicY<TF1>(g1->GetListOfFunctions());
720 
721  for (Int_t i = 0; i != g1->GetN(); ++i) {
722  g1->GetY()[i] = pow(10.0, g1->GetY()[i]);
723  }
724  }
725  }
726 
727 
728  /**
729  * Make x-axis and associated functions of given graph logarithmic (e.g.\ after using <tt>log10()</tt>).
730  *
731  * \param g1 graph
732  */
733  inline void setLogarithmicX(TGraphErrors* g1)
734  {
735  if (g1 != NULL) {
736 
737  setLogarithmicX<TF1>(g1->GetListOfFunctions());
738 
739  for (Int_t i = 0; i != g1->GetN(); ++i) {
740  g1->GetEX()[i] = pow(10.0, g1->GetX()[i] + g1->GetEX()[i]) - pow(10.0, g1->GetX()[i]);
741  g1->GetX() [i] = pow(10.0, g1->GetX()[i]);
742  }
743  }
744  }
745 
746 
747  /**
748  * Make y-axis and associated functions of given graph logarithmic (e.g.\ after using <tt>log10()</tt>).
749  *
750  * \param g1 graph
751  */
752  inline void setLogarithmicY(TGraphErrors* g1)
753  {
754  if (g1 != NULL) {
755 
756  setLogarithmicY<TF1>(g1->GetListOfFunctions());
757 
758  for (Int_t i = 0; i != g1->GetN(); ++i) {
759  g1->GetEY()[i] = pow(10.0, g1->GetY()[i] + g1->GetEY()[i]) - pow(10.0, g1->GetY()[i]);
760  g1->GetY() [i] = pow(10.0, g1->GetY()[i]);
761  }
762  }
763  }
764 
765 
766  /**
767  * Make x-axis and associated functions of given graph logarithmic (e.g.\ after using <tt>log10()</tt>).
768  *
769  * \param g2 graph
770  */
771  inline void setLogarithmicX(TGraph2D* g2)
772  {
773  if (g2 != NULL) {
774 
775  setLogarithmicX<TF2>(g2->GetListOfFunctions());
776 
777  for (Int_t i = 0; i != g2->GetN(); ++i) {
778  g2->GetX()[i] = pow(10.0, g2->GetX()[i]);
779  }
780  }
781  }
782 
783 
784  /**
785  * Make y-axis and associated functions of given graph logarithmic (e.g.\ after using <tt>log10()</tt>).
786  *
787  * \param g2 graph
788  */
789  inline void setLogarithmicY(TGraph2D* g2)
790  {
791  if (g2 != NULL) {
792 
793  setLogarithmicY<TF2>(g2->GetListOfFunctions());
794 
795  for (Int_t i = 0; i != g2->GetN(); ++i) {
796  g2->GetY()[i] = pow(10.0, g2->GetY()[i]);
797  }
798  }
799  }
800 
801 
802  /**
803  * Make x-axis and associated functions of given graph logarithmic (e.g.\ after using <tt>log10()</tt>).
804  *
805  * \param g2 graph
806  */
807  inline void setLogarithmicX(TGraph2DErrors* g2)
808  {
809  if (g2 != NULL) {
810 
811  setLogarithmicX<TF2>(g2->GetListOfFunctions());
812 
813  for (Int_t i = 0; i != g2->GetN(); ++i) {
814  g2->GetEX()[i] = pow(10.0, g2->GetX()[i] + g2->GetEX()[i]) - pow(10.0, g2->GetX()[i]);
815  g2->GetX() [i] = pow(10.0, g2->GetX()[i]);
816  }
817  }
818  }
819 
820 
821  /**
822  * Make y-axis and associated functions of given graph logarithmic (e.g.\ after using <tt>log10()</tt>).
823  *
824  * \param g2 graph
825  */
826  inline void setLogarithmicY(TGraph2DErrors* g2)
827  {
828  if (g2 != NULL) {
829 
830  setLogarithmicY<TF2>(g2->GetListOfFunctions());
831 
832  for (Int_t i = 0; i != g2->GetN(); ++i) {
833  g2->GetEY()[i] = pow(10.0, g2->GetY()[i] + g2->GetEY()[i]) - pow(10.0, g2->GetY()[i]);
834  g2->GetY() [i] = pow(10.0, g2->GetY()[i]);
835  }
836  }
837  }
838 
839 
840  /**
841  * Make x-axis of given multi-graph logarithmic (e.g.\ after using <tt>log10()</tt>).
842  *
843  * \param gn multi graph
844  */
845  inline void setLogarithmicX(TMultiGraph* gn)
846  {
847  if (gn != NULL) {
848  setLogarithmicX<TGraph>(gn->GetListOfGraphs());
849  }
850  }
851 
852 
853  /**
854  * Make y-axis of given multi-graph logarithmic (e.g.\ after using <tt>log10()</tt>).
855  *
856  * \param gn multi graph
857  */
858  inline void setLogarithmicY(TMultiGraph* gn)
859  {
860  if (gn != NULL) {
861  setLogarithmicY<TGraph>(gn->GetListOfGraphs());
862  }
863  }
864 
865 
866  /**
867  * Make x-axis of given line logarithmic (e.g.\ after using <tt>log10()</tt>).
868  *
869  * \param line line
870  */
871  inline void setLogarithmicX(TLine* line)
872  {
873  if (line != NULL) {
874  line->SetX1(pow(10.0, line->GetX1()));
875  line->SetX2(pow(10.0, line->GetX2()));
876  }
877  }
878 
879 
880  /**
881  * Make y-axis of given line logarithmic (e.g.\ after using <tt>log10()</tt>).
882  *
883  * \param line line
884  */
885  inline void setLogarithmicY(TLine* line)
886  {
887  if (line != NULL) {
888  line->SetY1(pow(10.0, line->GetY1()));
889  line->SetY2(pow(10.0, line->GetY2()));
890  }
891  }
892 
893 
894  /**
895  * Make x-axis of given ellipse logarithmic (e.g.\ after using <tt>log10()</tt>).
896  *
897  * \param ellipse ellipse
898  */
899  inline void setLogarithmicX(TEllipse* ellipse)
900  {
901  THROW(JFunctionalException, "Operation setLogarithmicX on TEllipse not allowed.");
902  }
903 
904 
905  /**
906  * Make y-axis of given ellipse logarithmic (e.g.\ after using <tt>log10()</tt>).
907  *
908  * \param ellipse ellipse
909  */
910  inline void setLogarithmicY(TEllipse* ellipse)
911  {
912  THROW(JFunctionalException, "Operation setLogarithmicY on TEllipse not allowed.");
913  }
914 
915 
916  /**
917  * Make x-axis of objects in list logarithmic (e.g.\ after using <tt>log10()</tt>).
918  *
919  * \param list list
920  */
921  template<class T>
922  inline void setLogarithmicX(TList* list)
923  {
924  for (TIter i(list); T* p = dynamic_cast<T*>(i.Next()); ) {
925  setLogarithmicX(p);
926  }
927  }
928 
929 
930  /**
931  * Make y-axis of objects in list logarithmic (e.g.\ after using <tt>log10()</tt>).
932  *
933  * \param list list
934  */
935  template<class T>
936  inline void setLogarithmicY(TList* list)
937  {
938  for (TIter i(list); T* p = dynamic_cast<T*>(i.Next()); ) {
939  setLogarithmicY(p);
940  }
941  }
942 
943 
944  /**
945  * Convert 1D histogram to PDF.
946  *
947  * Possible options are:
948  * - N normalise to histogram contents;
949  * - W divide by bin width;
950  * - E convert also bin errors.
951  *
952  * \param h1 histogram
953  * \param option option
954  * \param factor scaling factor
955  */
956  inline void convertToPDF(TH1& h1, const std::string& option = "NW", const double factor = 1.0)
957  {
958  using namespace std;
959 
960  const bool normalise = (option.find('N') != string::npos || option.find('n') != string::npos);
961  const bool bin_width = (option.find('W') != string::npos || option.find('w') != string::npos);
962  const bool use_error = (option.find('E') != string::npos || option.find('e') != string::npos);
963 
964  Double_t W = 1.0;
965 
966  if (normalise) {
967 
968  W = 0.0;
969 
970  for (Int_t i = 1; i <= h1.GetXaxis()->GetNbins(); ++i) {
971  W += h1.GetBinContent(i);
972  }
973  }
974 
975  if (W != 0.0) {
976 
977  for (Int_t i = 1; i <= h1.GetXaxis()->GetNbins(); ++i) {
978 
979  const Double_t w = W * (bin_width ? h1.GetXaxis()->GetBinWidth(i) : 1.0);
980 
981  h1.SetBinContent(i, h1.GetBinContent(i) * factor / w);
982 
983  if (use_error) {
984  h1.SetBinError(i, h1.GetBinError(i) * factor / w);
985  }
986  }
987  }
988  }
989 
990 
991  /**
992  * Convert 2D histogram to PDF.
993  *
994  * Possible options are:
995  * - N normalise to histogram contents;
996  * - X convert x-axis to PDF;
997  * - Y convert y-axis to PDF;
998  * - W divide by bin width;
999  * - E convert also bin errors.
1000  *
1001  * \param h2 histogram
1002  * \param option option
1003  * \param factor scaling factor
1004  */
1005  inline void convertToPDF(TH2& h2, const std::string& option = "NXYW", const double factor = 1.0)
1006  {
1007  using namespace std;
1008 
1009  const bool normalise = (option.find('N') != string::npos || option.find('n') != string::npos);
1010  const bool X = (option.find('X') != string::npos || option.find('x') != string::npos);
1011  const bool Y = (option.find('Y') != string::npos || option.find('y') != string::npos);
1012  const bool bin_width = (option.find('W') != string::npos || option.find('w') != string::npos);
1013  const bool use_error = (option.find('E') != string::npos || option.find('e') != string::npos);
1014 
1015  Double_t W = 1.0;
1016 
1017  if (X && Y) {
1018 
1019  if (normalise) {
1020 
1021  W = 0.0;
1022 
1023  for (Int_t ix = 1; ix <= h2.GetXaxis()->GetNbins(); ++ix) {
1024  for (Int_t iy = 1; iy <= h2.GetYaxis()->GetNbins(); ++iy) {
1025  W += h2.GetBinContent(ix,iy);
1026  }
1027  }
1028  }
1029 
1030  if (W != 0.0) {
1031 
1032  for (Int_t ix = 1; ix <= h2.GetXaxis()->GetNbins(); ++ix) {
1033  for (Int_t iy = 1; iy <= h2.GetYaxis()->GetNbins(); ++iy) {
1034 
1035  const Double_t w = W * (bin_width ? h2.GetXaxis()->GetBinWidth(ix) * h2.GetYaxis()->GetBinWidth(iy) : 1.0);
1036 
1037  h2.SetBinContent(ix, iy, h2.GetBinContent(ix,iy) * factor / w);
1038 
1039  if (use_error) {
1040  h2.SetBinError(ix, iy, h2.GetBinError(ix,iy) * factor / w);
1041  }
1042  }
1043  }
1044  }
1045 
1046  } else if (X) {
1047 
1048  for (Int_t iy = 1; iy <= h2.GetYaxis()->GetNbins(); ++iy) {
1049 
1050  if (normalise) {
1051 
1052  W = 0.0;
1053 
1054  for (Int_t ix = 1; ix <= h2.GetXaxis()->GetNbins(); ++ix) {
1055  W += h2.GetBinContent(ix,iy);
1056  }
1057  }
1058 
1059  if (W != 0.0) {
1060 
1061  for (Int_t ix = 1; ix <= h2.GetXaxis()->GetNbins(); ++ix) {
1062 
1063  const Double_t w = W * (bin_width ? h2.GetXaxis()->GetBinWidth(ix) : 1.0);
1064 
1065  h2.SetBinContent(ix, iy, h2.GetBinContent(ix,iy) * factor / w);
1066 
1067  if (use_error) {
1068  h2.SetBinError(ix, iy, h2.GetBinError(ix,iy) * factor / w);
1069  }
1070  }
1071  }
1072  }
1073 
1074  } else if (Y) {
1075 
1076  for (Int_t ix = 1; ix <= h2.GetXaxis()->GetNbins(); ++ix) {
1077 
1078  if (normalise) {
1079 
1080  W = 0.0;
1081 
1082  for (Int_t iy = 1; iy <= h2.GetYaxis()->GetNbins(); ++iy) {
1083  W += h2.GetBinContent(ix,iy);
1084  }
1085  }
1086 
1087  if (W != 0.0) {
1088 
1089  for (Int_t iy = 1; iy <= h2.GetYaxis()->GetNbins(); ++iy) {
1090 
1091  const Double_t w = W * (bin_width ? h2.GetYaxis()->GetBinWidth(iy) : 1.0);
1092 
1093  h2.SetBinContent(ix, iy, h2.GetBinContent(ix,iy) / w);
1094 
1095  if (use_error) {
1096  h2.SetBinError(ix, iy, h2.GetBinError(ix,iy) / w);
1097  }
1098  }
1099  }
1100  }
1101  }
1102  }
1103 
1104 
1105  /**
1106  * Convert 3D histogram to PDF.
1107  *
1108  * Possible options are:
1109  * - N normalise to histogram contents;
1110  * - X convert x-axis to PDF;
1111  * - Y convert y-axis to PDF;
1112  * - Z convert z-axis to PDF;
1113  * - W divide by bin width;
1114  * - E convert also bin errors.
1115  *
1116  * \param h3 histogram
1117  * \param option option
1118  * \param factor scaling factor
1119  */
1120  inline void convertToPDF(TH3& h3, const std::string& option = "NXYW", const double factor = 1.0)
1121  {
1122  using namespace std;
1123 
1124  const bool normalise = (option.find('N') != string::npos || option.find('n') != string::npos);
1125  const bool X = (option.find('X') != string::npos || option.find('x') != string::npos);
1126  const bool Y = (option.find('Y') != string::npos || option.find('y') != string::npos);
1127  const bool Z = (option.find('Z') != string::npos || option.find('z') != string::npos);
1128  const bool bin_width = (option.find('W') != string::npos || option.find('w') != string::npos);
1129  const bool use_error = (option.find('E') != string::npos || option.find('e') != string::npos);
1130 
1131  Double_t W = 1.0;
1132 
1133  if (X && Y && Z) {
1134 
1135  if (normalise) {
1136 
1137  W = 0.0;
1138 
1139  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1140  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1141  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1142  W += h3.GetBinContent(ix,iy,iz);
1143  }
1144  }
1145  }
1146  }
1147 
1148  if (W != 0.0) {
1149 
1150  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1151  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1152  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1153 
1154  const Double_t w = W * (bin_width ? h3.GetXaxis()->GetBinWidth(ix) * h3.GetYaxis()->GetBinWidth(iy) * h3.GetZaxis()->GetBinWidth(iz) : 1.0);
1155 
1156  h3.SetBinContent(ix, iy, iz, h3.GetBinContent(ix,iy,iz) * factor / w);
1157 
1158  if (use_error) {
1159  h3.SetBinError(ix, iy, iz, h3.GetBinError(ix,iy,iz) * factor / w);
1160  }
1161  }
1162  }
1163  }
1164  }
1165 
1166  } else if (X && Z) {
1167 
1168  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1169 
1170  if (normalise) {
1171 
1172  W = 0.0;
1173 
1174  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1175  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1176  W += h3.GetBinContent(ix,iy,iz);
1177  }
1178  }
1179  }
1180 
1181  if (W != 0.0) {
1182 
1183  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1184  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1185 
1186  const Double_t w = W * (bin_width ? h3.GetXaxis()->GetBinWidth(ix) * h3.GetZaxis()->GetBinWidth(iz) : 1.0);
1187 
1188  h3.SetBinContent(ix, iy, iz, h3.GetBinContent(ix,iy,iz) * factor / w);
1189 
1190  if (use_error) {
1191  h3.SetBinError(ix, iy, iz, h3.GetBinError(ix,iy,iz) * factor / w);
1192  }
1193  }
1194  }
1195  }
1196  }
1197 
1198  } else if (Y && Z) {
1199 
1200  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1201 
1202  if (normalise) {
1203 
1204  W = 0.0;
1205 
1206  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1207  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1208  W += h3.GetBinContent(ix,iy,iz);
1209  }
1210  }
1211  }
1212 
1213  if (W != 0.0) {
1214 
1215  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1216  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1217 
1218  const Double_t w = W * (bin_width ? h3.GetYaxis()->GetBinWidth(iy) * h3.GetZaxis()->GetBinWidth(iz) : 1.0);
1219 
1220  h3.SetBinContent(ix, iy, iz, h3.GetBinContent(ix,iy,iz) * factor / w);
1221 
1222  if (use_error) {
1223  h3.SetBinError(ix, iy, iz, h3.GetBinError(ix,iy,iz) * factor / w);
1224  }
1225  }
1226  }
1227  }
1228  }
1229 
1230  } else if (X && Y) {
1231 
1232  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1233 
1234  if (normalise) {
1235 
1236  W = 0.0;
1237 
1238  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1239  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1240  W += h3.GetBinContent(ix,iy,iz);
1241  }
1242  }
1243  }
1244 
1245  if (W != 0.0) {
1246 
1247  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1248  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1249 
1250  const Double_t w = W * (bin_width ? h3.GetXaxis()->GetBinWidth(ix) * h3.GetYaxis()->GetBinWidth(iy) : 1.0);
1251 
1252  h3.SetBinContent(ix, iy, iz, h3.GetBinContent(ix,iy,iz) * factor / w);
1253 
1254  if (use_error) {
1255  h3.SetBinError(ix, iy, iz, h3.GetBinError(ix,iy,iz) * factor / w);
1256  }
1257  }
1258  }
1259  }
1260  }
1261 
1262  } else if (X) {
1263 
1264  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1265  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1266 
1267  if (normalise) {
1268 
1269  W = 0.0;
1270 
1271  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1272  W += h3.GetBinContent(ix,iy,iz);
1273  }
1274  }
1275 
1276  if (W != 0.0) {
1277 
1278  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1279 
1280  const Double_t w = W * (bin_width ? h3.GetXaxis()->GetBinWidth(ix) : 1.0);
1281 
1282  h3.SetBinContent(ix, iy, iz, h3.GetBinContent(ix,iy,iz) * factor / w);
1283 
1284  if (use_error) {
1285  h3.SetBinError(ix, iy, iz, h3.GetBinError(ix,iy,iz) * factor / w);
1286  }
1287  }
1288  }
1289  }
1290  }
1291 
1292  } else if (Y) {
1293 
1294  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1295  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1296 
1297  if (normalise) {
1298 
1299  W = 0.0;
1300 
1301  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1302  W += h3.GetBinContent(ix,iy,iz);
1303  }
1304  }
1305 
1306  if (W != 0.0) {
1307 
1308  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1309 
1310  const Double_t w = W * (bin_width ? h3.GetYaxis()->GetBinWidth(iy) : 1.0);
1311 
1312  h3.SetBinContent(ix, iy, iz, h3.GetBinContent(ix,iy,iz) / w);
1313 
1314  if (use_error) {
1315  h3.SetBinError(ix, iy, iz, h3.GetBinError(ix,iy,iz) / w);
1316  }
1317  }
1318  }
1319  }
1320  }
1321 
1322  } else if (Z) {
1323 
1324  for (Int_t ix = 1; ix <= h3.GetXaxis()->GetNbins(); ++ix) {
1325  for (Int_t iy = 1; iy <= h3.GetYaxis()->GetNbins(); ++iy) {
1326 
1327  if (normalise) {
1328 
1329  W = 0.0;
1330 
1331  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1332  W += h3.GetBinContent(ix,iy,iz);
1333  }
1334  }
1335 
1336  if (W != 0.0) {
1337 
1338  for (Int_t iz = 1; iz <= h3.GetZaxis()->GetNbins(); ++iz) {
1339 
1340  const Double_t w = W * (bin_width ? h3.GetZaxis()->GetBinWidth(iz) : 1.0);
1341 
1342  h3.SetBinContent(ix, iy, iz, h3.GetBinContent(ix,iy,iz) / w);
1343 
1344  if (use_error) {
1345  h3.SetBinError(ix, iy, iz, h3.GetBinError(ix,iy,iz) / w);
1346  }
1347  }
1348  }
1349  }
1350  }
1351  }
1352  }
1353 
1354 
1355  /**
1356  * Set limits of TGraph.
1357  *
1358  * \param g1 graph
1359  */
1360  inline void setLimits(TGraph& g1)
1361  {
1362  using namespace std;
1363 
1364  Double_t ymin = numeric_limits<Double_t>::max();
1365  Double_t ymax = numeric_limits<Double_t>::lowest();
1366 
1367  for (Int_t i = 0; i != g1.GetN(); ++i) {
1368 
1369  const Double_t y = g1.GetY()[i];
1370 
1371  if (y > ymax) { ymax = y; }
1372  if (y < ymin) { ymin = y; }
1373  }
1374 
1375  g1.SetMinimum(ymin);
1376  g1.SetMaximum(ymax);
1377  }
1378 
1379 
1380  /**
1381  * Set limits of TGraph2D.
1382  *
1383  * \param g2 graph
1384  */
1385  inline void setLimits(TGraph2D& g2)
1386  {
1387  using namespace std;
1388 
1389  Double_t zmin = numeric_limits<Double_t>::max();
1390  Double_t zmax = numeric_limits<Double_t>::lowest();
1391 
1392  for (Int_t i = 0; i != g2.GetN(); ++i) {
1393 
1394  const Double_t z = g2.GetZ()[i];
1395 
1396  if (z > zmax) { zmax = z; }
1397  if (z < zmin) { zmin = z; }
1398  }
1399 
1400  g2.SetMinimum(zmin);
1401  g2.SetMaximum(zmax);
1402  }
1403 
1404 
1405  /**
1406  * Set axis range.
1407  *
1408  * \param xmin lower limit (I/O)
1409  * \param xmax upper limit (I/O)
1410  * \param logx logarithmic
1411  */
1412  inline void setRange(double& xmin,
1413  double& xmax,
1414  const bool logx)
1415  {
1416  if (logx) {
1417  xmin = log(xmin);
1418  xmax = log(xmax);
1419  }
1420 
1421  double dx = (xmax - xmin) * 0.1;
1422 
1423  if (logx || xmin >= dx || xmin < 0.0)
1424  xmin -= dx;
1425  else
1426  xmin = 0.0;
1427 
1428  xmax += dx;
1429 
1430  if (logx) {
1431  xmin = exp(xmin);
1432  xmax = exp(xmax);
1433  }
1434  }
1435 
1436 
1437  /**
1438  * Set axis with PMT address labels.
1439  *
1440  * This methods sets the labels of the given axis to the sorted values of the PMT ring and position.\n
1441  * It should normally be called before filling of the corresponding histogram.\n
1442  * The filling should then be made with the textual representation of the PMT ring and position (i.e.\ JDETECTOR::JPMTPhysicalAddress::toString).
1443  *
1444  * Alternatively, the filling can be made with the index of the PMT in the address map of the corresponding module
1445  * (i.e.\ JDETECTOR::JModuleAddressMap::getIndex).\n
1446  * In that case, the method can be called before or after filling of the histogram.
1447  *
1448  * \param axis axis
1449  * \param memo module address map
1450  */
1451  inline void setAxisLabels(TAxis *axis, const JModuleAddressMap& memo)
1452  {
1453  using namespace JPP;
1454 
1455  if (axis->GetNbins() == (int) memo.size()) {
1456 
1457  for (int i = 0; i != axis->GetNbins(); ++i) {
1458 
1459  const JPMTPhysicalAddress& address = memo[i];
1460 
1461  axis->SetBinLabel(i + 1, address.toString().c_str());
1462  }
1463 
1464  } else {
1465 
1466  THROW(JValueOutOfRange, "Number of bins " << axis->GetNbins() << " != " << memo.size());
1467  }
1468  }
1469 
1470 
1471  /**
1472  * Set axis labels with PMT addresses.
1473  *
1474  * This methods sets the labels of the given axis to the sorted values of the PMT ring and position.\n
1475  * It should be called after filling of the corresponding histogram.\n
1476  * The filling should have been made with the PMT number (e.g.\ KM3NETDAQ::JDAQHit::getPMT).
1477  *
1478  * \param h1 histogram
1479  * \param axis axis <tt>(x, X, y, Y, z, Z)</tt>
1480  * \param memo module address map
1481  */
1482  inline void setAxisLabels(TH1& h1, const std::string axis, const JModuleAddressMap& memo)
1483  {
1484  using namespace JPP;
1485 
1486  TAxis* pax = NULL;
1487 
1488  if (axis == "x" || axis == "X")
1489  pax = h1.GetXaxis();
1490  else if (axis == "y" || axis == "Y")
1491  pax = h1.GetYaxis();
1492  else if (axis == "z" || axis == "Z")
1493  pax = h1.GetZaxis();
1494  else
1495  THROW(JValueOutOfRange, "Invalid axis " << axis);
1496 
1497  if (pax->GetNbins() == (int) memo.size()) {
1498 
1499  for (int i = 0; i != pax->GetNbins(); ++i) {
1500 
1501  const JPMTPhysicalAddress& address = memo.getAddressTranslator(i);
1502 
1503  pax->SetBinLabel(i + 1, address.toString().c_str());
1504  }
1505 
1506  h1.LabelsOption("a", axis.c_str()); // sort labels
1507 
1508  } else {
1509 
1510  THROW(JValueOutOfRange, "Number of bins " << pax->GetNbins() << " != " << memo.size());
1511  }
1512  }
1513 
1514 
1515  /**
1516  * Check if given key corresponds to a TObject.
1517  *
1518  * \param key ROOT key
1519  * \return true if given key corresponds to a TObject; else false
1520  */
1521  inline bool isTObject(const TKey* key)
1522  {
1523  return (key != NULL && TClass::GetClass(key->GetClassName())->IsTObject());
1524  }
1525 
1526 
1527  /**
1528  * Get first object of which name matches given reguar expression.
1529  *
1530  * \param ls pointer to list of objects
1531  * \param name regular expression
1532  * \return pointer to object (maybe NULL)
1533  */
1534  inline TObject* getObject(TList* ls, const char* const name)
1535  {
1536  const TRegexp regexp(name);
1537 
1538  TIter iter(ls);
1539 
1540  for (TObject* p1; (p1 = (TObject*) iter.Next()) != NULL; ) {
1541 
1542  const TString tag(p1->GetName());
1543 
1544  if (tag.Contains(regexp)) {
1545  return p1;
1546  }
1547  }
1548 
1549  return NULL;
1550  }
1551 
1552 
1553  /**
1554  * Get function.
1555  *
1556  * \param h1 histogram
1557  * \param fcn function name
1558  * \return pointer to function (maybe NULL)
1559  */
1560  inline TF1* getFunction(TH1* h1, const char* const fcn)
1561  {
1562  return (TF1*) getObject(h1->GetListOfFunctions(), fcn);
1563  }
1564 
1565 
1566  /**
1567  * Get function.
1568  *
1569  * \param g1 graph
1570  * \param fcn function name
1571  * \return pointer to function (maybe NULL)
1572  */
1573  inline TF1* getFunction(TGraph* g1, const char* const fcn)
1574  {
1575  return (TF1*) getObject(g1->GetListOfFunctions(), fcn);
1576  }
1577 
1578 
1579  /**
1580  * Get function.
1581  *
1582  * \param g2 graph
1583  * \param fcn function name
1584  * \return pointer to function (maybe NULL)
1585  */
1586  inline TF1* getFunction(TGraph2D* g2, const char* const fcn)
1587  {
1588  return (TF1*) getObject(g2->GetListOfFunctions(), fcn);
1589  }
1590 }
1591 
1592 #endif
static const char *const sqrt()
Set contents to signed difference between squares.
const double xmax
Definition: JQuadrature.cc:24
data_type w[N+1][M+1]
Definition: JPolint.hh:867
Exceptions.
TObject * getObject(const JRootObjectID &id)
Get first TObject with given identifier.
double getValue(const JScale_t scale)
Get numerical value corresponding to scale.
Definition: JScale.hh:47
TF1 * getFunction(TH1 *h1, const char *const fcn)
Get function.
int getParameter(const std::string &text)
Get parameter number from text string.
static const char *const Divide()
ROOT TH1::Divide.
TPaveText * p1
static const char *const Subtract()
ROOT TH1::Subtract.
Exception for a functional operation.
Definition: JException.hh:142
static const char *const SAME_AS_INPUT()
Set name of output histogram to name of input histogram.
then usage $script[< detector identifier >< run range >]< QA/QCfile > nExample script to produce data quality plots nWhen a detector identifier and run range are data are downloaded from the database nand subsequently stored in the given QA QC file
Definition: JDataQuality.sh:19
char text[TEXT_SIZE]
Definition: elog.cc:72
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:712
Definition: JRoot.hh:19
static const char *const Replace()
Set contents to associated function.
void setLimits(TGraph &g1)
Set limits of TGraph.
static const char *const subtract()
Subtract contents with lookup bin in second histogram.
Auxiliary class to handle file name, ROOT directory and object name.
then fatal Wrong number of arguments fi set_variable STRING $argv[1] set_variable DETECTORXY_TXT $WORKDIR $DETECTORXY_TXT tail read X Y CHI2 RMS printf optimum n $X $Y $CHI2 $RMS awk v Y
static const char *const stdev()
Set contents to standard deviation.
is
Definition: JDAQCHSM.chsm:167
then echo The file $DIR KM3NeT_00000001_00000000 root already please rename or remove it first
void setLogarithmic(TAxis *axis)
Make histogram axis logarithmic (e.g. after using log10()).
const JPolynome f1(1.0, 2.0, 3.0)
Function.
bool isTAttLine(const TObject *object)
Get drawing option of TH1.
static const char *const efficiency()
Divide contents and multiply errors with inefficiency.
static const char *const Add()
ROOT TH1::Add.
Lookup table for PMT addresses in optical module.
static const char *const multiply()
Multiply contents with lookup bin in second histogram.
TFile * getFile(const std::string &file_name, const std::string &option="exist")
Get TFile pointer corresponding to give file name.
static const char *const add()
Add contents with lookup bin in second histogram.
Auxiliary data structure for JOpera1D.cc and JOpera2D.cc applications.
void setLogarithmicX(TList *list)
Make x-axis of objects in list logarithmic (e.g. after using log10()).
Double_t getResult(const TString &text, TObject *object=NULL)
Get result of given textual formula.
do set_variable OUTPUT_DIRECTORY $WORKDIR T
void convertToPDF(TH1 &h1, const std::string &option="NW", const double factor=1.0)
Convert 1D histogram to PDF.
void setAxisLabels(TAxis *axis, const JModuleAddressMap &memo)
Set axis with PMT address labels.
T pow(const T &x, const double y)
Power .
Definition: JMath.hh:97
static const char *const Multiply()
ROOT TH1::Multiply.
static const char *const SAME_AS_OPERATION()
Set name of output histogram to name of operation.
then fatal The output file must have the wildcard in the name
Definition: JCanberra.sh:31
void setLogarithmicY(TList *list)
Make y-axis of objects in list logarithmic (e.g. after using log10()).
const double xmin
Definition: JQuadrature.cc:23
void setRange(double &xmin, double &xmax, const bool logx)
Set axis range.
const JPMTAddressTranslator & getAddressTranslator(const int tdc) const
Get PMT address translator.
then usage $script< input file >[option[primary[working directory]]] nWhere option can be N
Definition: JMuonPostfit.sh:40
then fatal The output file must have the wildcard in the e g root fi eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY JAcoustics sh $DETECTOR_ID source JAcousticsToolkit sh CHECK_EXIT_CODE typeset A EMITTERS get_tripods $WORKDIR tripod txt EMITTERS get_transmitters $WORKDIR transmitter txt EMITTERS for EMITTER in
Definition: JCanberra.sh:48
static const char *const divide()
Divide contents with lookup bin in second histogram.
void copy(const Head &from, JHead &to)
Copy header from from to to.
Definition: JHead.cc:162
Exception for parsing value.
Definition: JException.hh:196
bool isTObject(const TKey *key)
Check if given key corresponds to a TObject.
no fit printf nominal n $STRING awk v X
TDirectory * getDirectory(const JRootObjectID &id)
Get TDirectory pointer.
std::string getFilename(const std::string &file_name)
Get file name part, i.e. part after last JEEP::PATHNAME_SEPARATOR if any.
Definition: JeepToolkit.hh:128
do set_variable MODULE getModule a $WORKDIR detector_a datx L $STRING JEditDetector a $WORKDIR detector_a datx M $MODULE setz o $WORKDIR detector_a datx JEditDetector a $WORKDIR detector_b datx M $MODULE setz o $WORKDIR detector_b datx done echo Output stored at $WORKDIR detector_a datx and $WORKDIR tripod_a txt JDrawDetector2D a $WORKDIR detector_a datx a $WORKDIR detector_b datx L BL o detector $FORMAT $BATCH JDrawDetector2D T $WORKDIR tripod_a txt T $WORKDIR tripod_b txt L BL o tripod $FORMAT $BATCH JCompareDetector a $WORKDIR detector_a datx b $WORKDIR detector_b datx o $WORKDIR abc root &dev null for KEY in X Y Z
Exception for accessing a value in a collection that is outside of its range.
Definition: JException.hh:178
TString getLogarithmic(const TString &formula, const char parameter)
Make given parameter in formula logarithmic (e.g. after using log10()).
static const char *const TIMESTAMP
Time stamp of earliest UTC time.
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable STRING $argv[2] set_array QUANTILES set_variable FORMULA *[0] exp(-0.5 *(x-[1])*(x-[1])/([2]*[2]))" set_variable MODULE `getModule -a $DETECTOR -L "$STRING 0"` source JAcousticsToolkit.sh typeset -A TRIPODS get_tripods $WORKDIR/tripod.txt TRIPODS XMEAN
Double_t g1(const Double_t x)
Function.
Definition: JQuantiles.cc:25