19 This module contains set of simple and useful utilities for booking and
20 manipulations with Gaudi-AIDA histograms, inspired by Thomas' request
22 The module contains following public symbols:
24 - book for booking of various 1D,2D&3D-histograms
25 - bookProf for booking of various 1D&2D-profiles
26 - getAsAIDA for retrieval of histograms/profiles from HTS in AIDA format
27 - getAsROOT for retrieval of histograms/profiles from HTS in ROOT format
28 - fill for smart filling of 1D-histogram (AIDA or ROOT)
29 - aida2root for conversion of AIDA histogram to ROOT
30 - HistoStats for many statistical information
31 - HistoFile class for storing histograms to a file
32 - histoDump for dumping of the histogram in text format (a'la HBOOK)
33 - dumpHisto for dumping of the histogram in text format (a'la HBOOK)
36 from __future__
import print_function
39 __author__ =
"Vanya BELYAEV ibelyaev@physics.syr.edu"
71 Helper private auxiliary function to get Application Manager
73 gaudi = kwargs.get(
"gaudi",
None)
77 raise RuntimeError(
"Unable to get valid ApplicationMgr")
79 state = gaudi._isvc.FSMState()
80 if state < cpp.Gaudi.StateMachine.CONFIGURED:
82 state = gaudi._isvc.FSMState()
83 if state < cpp.Gaudi.StateMachine.INITIALIZED:
95 Helper private auxiliary function to get iHistogramSvs
97 svc = kwargs.get(
"service",
None)
99 svc = kwargs.get(
"svc",
None)
103 return gaudi.histsvc()
112 Helper private auxiliary function to get iDataSvs
114 svc = kwargs.get(
"service",
None)
116 svc = kwargs.get(
"svc",
None)
120 return gaudi.evtsvc()
129 The trivial function to book the various 1D,2D&3D-histograms
131 (1) book the trivial 1D histogram with full path
133 >>> h1D = book ( 'path/to/my/histo' , ## path in Histogram Transient Store
134 'cosine of decay angle ' , ## histogram title
135 100 , ## number of bins
139 (2) book the trivial 1D histogram with directory path and string ID :
141 >>> h1D = book ( 'path/to/directory' , ## the path to directory in HTS
142 'H1' , ## string histogram identifier
143 'cosine of decay angle ' , ## histogram title
144 100 , ## number of bins
148 (3) book the trivial 1D histogram with directory path and integer ID :
150 >>> h1D = book ( 'path/to/directory' , ## the path to directory in HTS
151 124 , ## integer histogram identifier
152 'cosine of decay angle ' , ## histogram title
153 100 , ## number of bins
157 (4) book the trivial 2D histogram with full path
159 >>> h1D = book ( 'path/to/my/histo' , ## path in Histogram Transient Store
160 'm12**2 versus m13**2' , ## histogram title
161 50 , ## number of X-bins
164 50 , ## number of Y-bins
168 (5) book the trivial 2D histogram with directory path and literal ID
170 >>> h1D = book ( 'path/to/directory' , ## path in Histogram Transient Store
171 'Dalitz1' , ## literal histogram identifier
172 'm12**2 versus m13**2' , ## histogram title
173 50 , ## number of X-bins
176 50 , ## number of Y-bins
180 (6) book the trivial 2D histogram with directory path and integer ID
182 >>> h1D = book ( 'path/to/directory' , ## path in Histogram Transient Store
183 854 , ## integer histogram identifier
184 'm12**2 versus m13**2' , ## histogram title
185 50 , ## number of X-bins
188 50 , ## number of Y-bins
192 (7) book the trivial 3D histogram with full path
194 >>> h1D = book ( 'path/to/my/histo' , ## path in Histogram Transient Store
195 'x vs y vs z ' , ## histogram title
196 10 , ## number of X-bins
199 10 , ## number of Y-bins
202 10 , ## number of Z-bins
206 (8) book the trivial 3D histogram with directory path and literal ID
208 >>> h1D = book ( 'path/to/directory' , ## path in Histogram Transient Store
209 'xyz' , ## literal histogram identifier
210 'x vs y vs z' , ## histogram title
211 10 , ## number of X-bins
214 10 , ## number of Y-bins
217 10 , ## number of Z-bins
221 (9) book the trivial 3D histogram with directory path and integer ID
223 >>> h1D = book ( 'path/to/directory' , ## path in Histogram Transient Store
224 888 , ## integer histogram identifier
225 'x vs y vs z' , ## histogram title
226 10 , ## number of X-bins
229 10 , ## number of Y-bins
232 10 , ## number of Z-bins
236 Many other booking methods are available,
237 e.g. for the histograms with non-equidistant bins, see IHistogamSvc::book
240 if useROOT
or kwargs.get(
"useROOT",
False)
or not kwargs.get(
"useAIDA",
True):
241 from ROOT
import TH1D
246 if not str
is type(a1):
249 return TH1D(a0 + a1, a2, *args[3:])
251 return TH1D(a0, a1, *args[2:])
255 raise RuntimeError(
"Unable to get valid HistogramService ")
257 return svc.book(*args)
261 "\n\n" +
"\thelp(iHistogramSvc.book) : \n\n" + iHistogramSvc.book.__doc__
264 "\n\n" +
"\thelp(IHistogramSvc::book) : \n\n" + cpp.IHistogramSvc.book.__doc__
274 The trivial function to book 1D&2D profile histograms:
276 (1) book 1D-profile histogram with full path in Histogram Transient Store:
277 >>> histo = bookProf ( 'path/to/my/profile' , ## path in Histogram Transient Store
278 'Energy Correction' , ## the histogram title
279 100 , ## number of X-bins
283 (2) book 1D-profile histogram with directory path and literal ID
284 >>> histo = bookProf ( 'path/to/my/profile' , ## path in Histogram Transient Store
285 'Calibration' , ## the histogram literal identifier
286 'Energy Correction' , ## the histogram title
287 100 , ## number of X-bins
291 (3) book 1D-profile histogram with directory path and integer ID
292 >>> histo = bookProf ( 'path/to/my/profile' , ## path in Histogram Transient Store
293 418 , ## the histogram integer identifier
294 'Energy Correction' , ## the histogram title
295 100 , ## number of X-bins
299 (4) book 2D-profile histogram with full path in Histogram Transient Store:
300 >>> histo = bookProf ( 'path/to/my/profile' , ## path in Histogram Transient Store
301 'Energy Correction' , ## the histogram title
302 50 , ## number of X-bins
305 50 , ## number of Y-bins
309 (5) book 2D-profile histogram with directory path and literal ID
310 >>> histo = bookProf ( 'path/to/my/profile' , ## path in Histogram Transient Store
311 'Calibration' , ## the histogram literal identifier
312 'Energy Correction' , ## the histogram title
313 50 , ## number of X-bins
316 50 , ## number of Y-bins
320 (6) book 2D-profile histogram with directory path and integer ID
321 >>> histo = bookProf ( 'path/to/my/profile' , ## path in Histogram Transient Store
322 418 , ## the histogram integer identifier
323 'Energy Correction' , ## the histogram title
324 50 , ## number of X-bins
327 50 , ## number of Y-bins
331 Many other booking methods are available,
332 e.g. for the histograms with non-equidistant bins, see IHistogamSvc::book
337 raise RuntimeError(
"Unable to get valid HistogramService ")
339 return svc.bookProf(*args)
342 bookProf.__doc__ += (
343 "\n\n" +
"\thelp(iHistogramSvc.bookProf) : \n\n" + iHistogramSvc.bookProf.__doc__
345 bookProf.__doc__ += (
347 +
"\thelp(IHistogramSvc::bookProf) : \n\n"
348 + cpp.IHistogramSvc.bookProf.__doc__
358 The most trivial function to retrieve the histogram from Histogram Transient Store
359 The histogram is returned by reference to its AIDA-representation (if possible)
361 >>> h = getAsAIDA ( 'some/path/to/my/histogram' )
366 raise RuntimeError(
"Unable to get valid HistogramService ")
368 return svc.getAsAIDA(path)
371 getAsAIDA.__doc__ += (
372 "\n\n" +
"\thelp(iHistogramSvc.getAsAIDA) : \n\n" + iHistogramSvc.getAsAIDA.__doc__
374 getAsAIDA.__doc__ += (
375 "\n\n" +
"\thelp(iHistogramSvc.retrieve) : \n\n" + iHistogramSvc.retrieve.__doc__
385 The most trivial function to retrieve the histogram from Histogram Transient Store
386 The histogram is returned by reference to its underlying native ROOT-representation (if possible)
388 >>> h = getAsROOT ( 'some/path/to/my/histogram' )
393 raise RuntimeError(
"Unable to get valid HistogramService ")
395 return svc.getAsROOT(path)
398 getAsROOT.__doc__ += (
399 "\n\n" +
"\thelp(iHistogramSvc.getAsROOT) : \n\n" + iHistogramSvc.getAsROOT.__doc__
414 The function which allows 'the smart fill' of 1D-histogram
418 The most trivial case, fill with the value
419 >>> fill ( histo , 1.0 )
421 Fill from any iterable object (sequence)
422 >>> fill ( histo , [1,,2,3,4,5,10] )
424 Fill from iterable object and apply the function:
425 >>> fill ( histo , [1,2,3,4,5] , math.sin )
428 >>> fill ( histo , [1,2,3,4,5] , lambda x : x*x )
431 >>> fill ( histo , [1,2,3,4,5] , fun = lambda x : x*x )
433 Use internal attributes:
434 >>> tracks = evtSvc['Rec/Track/Best'] ## iterable container of tracks
435 >>> fill ( histo , tracks , lambda t : t.pt() )
437 Apply the predicate: fill only even numbers:
438 >>> fill ( histo , [1,2,3,4,5,6,7] , lambda x : x , lambda y : y%2 )
440 The same (omit the trivial function) :
441 >>> fill ( histo , [1,2,3,4,5,6,7] , cut = lambda y : y%2 )
443 Apply the predicate: fill only pt for positively charged tracks:
444 >>> tracks = evtSvc['Rec/Track/Best']
445 >>> fill ( histo , tracks , lambda t : t.pt() , lambda t : 0<t.charge() )
448 >>> fill ( histo , tracks ,
449 fun = lambda t : t.pt() ,
450 cut = lambda t : 0<t.charge() )
452 Ordinary functions are also fine:
453 >>> def myfun ( track ) : return sin( track.pt() + track.p() )
454 >>> def mycut ( track ) : return track.p() > 100 * GeV
455 >>> fill ( histo , tracks , myfun , mycut )
457 The 'data' could be the address in TES, in this case the object
458 is retrieved from TES and the method is applied to the objects,
460 >>> fill ( histo , ## the reference to the histogram
461 'Rec/Track/Best' , ## the location of objects in TES
462 lambda t : t.pt() ) ## function to be used for histogram fill
463 >>> fill ( histo , ## the reference to the histogram
464 'Rec/Track/Best' , ## the address of objects in TES
465 lambda t : t.pt() , ## the function to be used for histogram fill
466 lambda t : t.charge()>0 ) ## the criteria to select tracks
468 The arguments 'fun' and 'cut' could be strings, in this case
469 they are evaluated by python before usage.
470 This option could be often very useful.
475 if type(data) == str:
478 return fill(histo, data, fun, cut, **kwargs)
482 fun = eval(fun, globals())
486 cut = eval(cut, globals())
488 if not hasattr(data,
"__iter__"):
491 if not hasattr(histo,
"fill")
and hasattr(histo,
"Fill"):
492 setattr(histo,
"fill", getattr(histo,
"Fill"))
497 histo.fill(fun(item))
504 aida2root = cpp.Gaudi.Utils.Aida2ROOT.aida2root
514 >>> aida = ... ## get AIDA histogram
515 >>> root = aida.toROOT() ## convert it to ROOT
521 _to_root_.__doc__ += aida2root.__doc__
524 cpp.AIDA.IHistogram3D,
525 cpp.AIDA.IHistogram2D,
526 cpp.AIDA.IHistogram1D,
530 if not hasattr(t,
"Fill")
and hasattr(t,
"fill"):
531 setattr(t,
"Fill", getattr(t,
"fill"))
532 for attr
in (
"toROOT",
"toRoot",
"asROOT",
"asRoot",
"AsROOT",
"AsRoot"):
533 if not hasattr(t, attr):
534 setattr(t, attr, _to_root_)
536 cpp.AIDA.IHistogram3D.__repr__ =
lambda s: cpp.GaudiAlg.Print3D.toString(
539 cpp.AIDA.IHistogram3D.__str__ = cpp.AIDA.IHistogram3D.__repr__
541 HistoStats = cpp.Gaudi.Utils.HistoStats
549 Evaluate 'bin-by-bin' momentum of order 'order' around the value 'value'
553 >>> print(h1.moment ( 5 ))
556 return HistoStats.moment(self, order, value)
565 Evaluate error for 'bin-by-bin' momentum of order 'order' around the value 'value'
569 >>> print(h1.momentErr ( 5 ))
572 return HistoStats.momentErr(self, order)
581 Evaluate 'bin-by-bin' central momentum (around mean value)
585 >>> print(h1.centralMoment ( 5 ))
588 return HistoStats.centralMoment(self, order)
597 Evaluate error for 'bin-by-bin' central momentum (around mean value)
601 >>> print(h1.centralMomentErr ( 5 ))
604 return HistoStats.centralMomentErr(self, order)
613 Evaluate 'bin-by-bin' skewness for 1D AIDA histogram
616 >>> print(h1.skewness())
619 return HistoStats.skewness(self)
628 Evaluate error for 'bin-by-bin' skewness
631 >>> print(h1.skewnessErr())
634 return HistoStats.skewnessErr(self)
643 Evaluate 'bin-by-bin' kurtosis
646 >>> print(h1.kurtosis ())
649 return HistoStats.kurtosis(self)
658 Evaluate error for 'bin-by-bin' kurtotis for 1D AIDA histogram
661 >>> print(h1.kurtotisErr())
664 return HistoStats.kurtosisErr(self)
672 Number of equivalent entries
674 return HistoStats.nEff(self)
682 Evaluate the MEAN value
684 return HistoStats.mean(self)
692 Evaluate the error for MEAN estimate
694 return HistoStats.meanErr(self)
702 Evaluate the RMS for AIDA histogram
704 return HistoStats.rms(self)
712 Evaluate the error for RMS estimate
714 return HistoStats.rmsErr(self)
722 Get an error in the sum bin height ('in-range integral')
724 return HistoStats.sumBinHeightErr(self)
731 """Get an error in the sum bin height ('in-range integral')"""
732 return HistoStats.sumAllBinHeightErr(self)
740 The fraction of overflow entries (useful for shape comparison)
742 return HistoStats.overflowEntriesFrac(self)
750 The error for fraction of overflow entries (useful for shape comparison)
752 return HistoStats.overflowEntriesFracErr(self)
760 The fraction of underflow entries (useful for shape comparison)
762 return HistoStats.underflowEntriesFrac(self)
770 The error for fraction of underflow entries (useful for shape comparison)
772 return HistoStats.underflowEntriesFracErr(self)
780 The fraction of overflow integral (useful for shape comparison)
782 return HistoStats.overflowIntegralFrac(self)
790 The error for fraction of overflow integral (useful for shape comparison)
792 return HistoStats.overflowIntegralFracErr(self)
800 The fraction of underflow integral (useful for shape comparison)
802 return HistoStats.underflowIntegralFrac(self)
810 The error for fraction of underflow integral (useful for shape comparison)
812 return HistoStats.underflowIntegralFracErr(self)
823 Get number of entries in histogram up to the certain bin (not-included)
825 attention: underflow bin is included!
828 >>> print(h1.nEntries ( 10 ))
830 Get number of entries in histogram form the certain
831 minimal bin up to the certain maximal bin (not-included)
834 >>> print(h1.nEntries ( 10 , 15 ))
837 if i2 < i1
or i2 < 0:
838 return HistoStats.nEntries(self, i1)
839 return HistoStats.nEntries(self, i1, i2)
847 Get the fraction of entries in histogram up to the certain bin (not-included)
849 attention: underflow bin is included!
852 >>> print(h1.nEntriesFrac ( 10 ))
854 Get the fraction of entries in histogram form the certain
855 minimal bin up to the certain maximal bin (not-included)
858 >>> print(h1.nEntriesFrac ( 10 , 15 ))
861 if i2 < i1
or i2 < 0:
862 return HistoStats.nEntriesFrac(self, i1)
863 return HistoStats.nEntriesFrac(self, i1, i2)
871 Get error for fraction of entries in histogram up to the certain bin (not-included)
873 attention: underflow bin is included!
876 >>> print(h1.nEntriesFracErr( 10 ))
878 Get error fraction of entries in histogram form the certain
879 minimal bin up to the certain maximal bin (not-included)
882 >>> print(h1.nEntriesFracErr ( 10 , 15 ))
885 if i2 < i1
or i2 < 0:
886 return HistoStats.nEntriesFracErr(self, i1)
887 return HistoStats.nEntriesFracErr(self, i1, i2)
891 i1DH = cpp.AIDA.IHistogram1D
893 if not hasattr(i1DH,
"moment"):
894 i1DH.moment = _moment_
895 if not hasattr(i1DH,
"momentErr"):
896 i1DH.momentErr = _momentErr_
897 if not hasattr(i1DH,
"centralMoment"):
898 i1DH.centralMoment = _centralMoment_
899 if not hasattr(i1DH,
"momentMomentErr"):
900 i1DH.centralMomentErr = _centralMomentErr_
901 if not hasattr(i1DH,
"nEff"):
903 if not hasattr(i1DH,
"mean"):
905 if not hasattr(i1DH,
"meanErr"):
906 i1DH.meanErr = _meanErr_
907 if not hasattr(i1DH,
"rms"):
909 if not hasattr(i1DH,
"rmsErr"):
910 i1DH.rmsErr = _rmsErr_
911 if not hasattr(i1DH,
"skewness"):
912 i1DH.skewness = _skewness_
913 if not hasattr(i1DH,
"skewnessErr"):
914 i1DH.skewnessErr = _skewnessErr_
915 if not hasattr(i1DH,
"kurtosis"):
916 i1DH.kurtosis = _kurtosis_
917 if not hasattr(i1DH,
"kurtosisErr"):
918 i1DH.kurtosisErr = _kurtosisErr_
920 if not hasattr(i1DH,
"overflowEntriesFrac"):
921 i1DH.overflowEntriesFrac = _overflowEntriesFrac_
922 if not hasattr(i1DH,
"overflowEntriesFracErr"):
923 i1DH.overflowEntriesFracErr = _overflowEntriesFracErr_
924 if not hasattr(i1DH,
"underflowEntriesFrac"):
925 i1DH.underflowEntriesFrac = _underflowEntriesFrac_
926 if not hasattr(i1DH,
"underflowEntriesFracErr"):
927 i1DH.underflowEntriesFracErr = _underflowEntriesFracErr_
929 if not hasattr(i1DH,
"overflowIntegralFrac"):
930 i1DH.overflowIntegralFrac = _overflowIntegralFrac_
931 if not hasattr(i1DH,
"overflowIntegralFracErr"):
932 i1DH.overflowIntegralFracErr = _overflowIntegralFracErr_
933 if not hasattr(i1DH,
"underflowIntegralFrac"):
934 i1DH.underflowIntegralFrac = _underflowIntegralFrac_
935 if not hasattr(i1DH,
"underflowIntegralFracErr"):
936 i1DH.underflowIntegralFracErr = _underflowIntegralFracErr_
938 if not hasattr(i1DH,
"nEntries"):
939 i1DH.nEntries = _nEntries_
940 if not hasattr(i1DH,
"nEntriesFrac"):
941 i1DH.nEntriesFrac = _nEntriesFrac_
942 if not hasattr(i1DH,
"nEntriesFracErr"):
943 i1DH.nEntriesFracErr = _nEntriesFracErr_
950 Get the path in THS for the given AIDA object:
953 >>> print(aida.path())
956 return cpp.Gaudi.Utils.Histos.path(self)
959 iBH = cpp.AIDA.IBaseHistogram
960 if not hasattr(iBH,
"path"):
962 if not hasattr(iBH,
"TESpath"):
964 if not hasattr(iBH,
"location"):
965 iBH.location = _path_
972 Dump the histogram/profile in text format (a'la HBOOK)
975 >>> print(dumpHisto ( histo ))
977 >>> print(histo.dump())
978 >>> print(histo.dump( 20 , 20 ))
979 >>> print(histo.dump( 20 , 20 , True ))
984 return cpp.Gaudi.Utils.Histos.histoDump(histo, *args)
987 __dumpHisto__.__doc__ =
"\n" + cpp.Gaudi.Utils.Histos.histoDump.__doc__
991 histoDump = __dumpHisto__
992 dumpHisto = __dumpHisto__
995 cpp.AIDA.IHistogram1D,
1002 for method
in (
"dump",
"dumpHisto",
"dumpAsText"):
1003 if not hasattr(t, method):
1004 setattr(t, method, __dumpHisto__)
1011 Class to write histograms to a ROOT file.
1012 hFile = HistoFile("myFile.root")
1018 hFile.save(myHisto0, myHisto1, myHisto2)
1019 histoList = [h0, h1, h2, h3]
1020 hFile.save(histoList)
1025 __author__ =
"Juan Palacios juan.palacios@nikhef.nl"
1028 self.
file = ROOT.TFile(fileName,
"RECREATE")
1029 from GaudiPython
import gbl
1033 gbl.AIDA.IHistogram1D,
1034 gbl.AIDA.IHistogram2D,
1035 gbl.AIDA.IHistogram3D,
1036 gbl.AIDA.IProfile1D,
1037 gbl.AIDA.IProfile2D,
1038 gbl.AIDA.IHistogram,
1042 histoType =
type(histo)
1050 This function stores histograms on the file for future saving.
1051 It takes an arbitrary number of AIDA or ROOT histograms or
1055 if type(args[0]) == list:
1070 if "__main__" == __name__:
1076 print(sys.modules[__name__].__dict__[o].__doc__)