4 """ Compare ROOT objects and perform e.g. chi2 tests.
5 A small command line interface for testing/debugging purposes is included.
6 Run `python3 validationcomparison.py --help` for more information. """
9 from abc
import ABC, abstractmethod
13 from typing
import Optional
19 from metaoptions
import MetaOptionParser
31 class ComparisonFailed(Exception):
33 The comparison failed for some reason. For example
34 because ROOT was not able to compute the Chi^2 properly
39 class ObjectsNotSupported(Exception):
41 The type and/or combination of provided ROOT objects
42 is not supported for comparison
47 class DifferingBinCount(Exception):
49 The two ROOT objects provided have a different bin count
50 and therefor, cannot be compared using the Chi2 test
55 class TooFewBins(Exception):
57 Not sufficient bins to perform the Chi^2 test
70 mop: Optional[MetaOptionParser]
71 ) ->
"ComparisonBase":
72 """ Uses the metaoptions to determine which comparison algorithm is used
73 and initializes the corresponding subclass of :class:`ComparisonBase` that
74 implements the actual comparison and holds the results.
75 @param object_1 ROOT TObject
76 @param object_2 ROOT TObject
77 @param mop Metaoption parser
79 if mop.has_option(
"kolmogorov"):
80 tester = KolmogorovTest
81 elif mop.has_option(
"andersondarling"):
82 tester = AndersonDarlingTest
99 class ComparisonBase(ABC):
101 Base class for all comparison implementations.
105 1. Initialize the class together with two ROOT objects of different
106 revisions (that are to be compared) and the metaoptions (given in the
107 corresponding validation (steering) file), that determine how to compare
110 2. The Comparison class saves the ROOT objects and the metaoptions
111 internally, but does not compute anything yet
113 3. If :meth:`ensure_compute` is called, or any property is accessed that
114 depends on computation, the internal implementation :meth:`_compute`
115 (to be implemented in the subclass) is called.
117 4. :meth:`_compute` ensures that all values, like chi2, p-value etc. are
120 5. Two properties :meth:`comparison_result` (pass/warning/error) and
121 :meth:`comparison_result_long` (longer description of the comparison result)
122 allow to access the results.
125 def __init__(self, object_a, object_b,
126 mop: Optional[MetaOptionParser] =
None, debug=
False):
128 Initialize ComparisonBase class
132 :param mop: MetaOptionParser
133 :param debug (bool): Debug mode enabled?
136 self.object_a = object_a
139 self.object_b = object_b
143 mop = MetaOptionParser()
150 self.computed =
False
153 self._comparison_result =
"not_compared"
156 self._comparison_result_long =
""
158 def ensure_compute(self):
160 Ensure all required quantities get computed and are cached inside the
166 if self.mop.has_option(
"nocompare"):
168 self._comparison_result_long =
'Testing is disabled for this plot'
171 fail_message =
"Comparison failed: "
176 except ObjectsNotSupported
as e:
177 self._comparison_result_long = fail_message + str(e)
178 except DifferingBinCount
as e:
179 self._comparison_result =
"error"
180 self._comparison_result_long = fail_message + str(e)
181 except TooFewBins
as e:
182 self._comparison_result_long = fail_message + str(e)
183 except ComparisonFailed
as e:
184 self._comparison_result =
"error"
185 self._comparison_result_long = fail_message + str(e)
186 except Exception
as e:
187 self._comparison_result =
"error"
188 self._comparison_result_long =
"Unknown error occurred. Please " \
189 "submit a bug report. " + str(e)
193 self._comparison_result_long = self._get_comparison_result_long()
194 self._comparison_result = self._get_comparison_result()
199 def _get_comparison_result(self) -> str:
200 """ Used to format the value of :attr:`_comparison_result`. """
204 def _get_comparison_result_long(self) -> str:
205 """ Used to format the value of :attr:`_comparison_result_long`. """
209 def comparison_result(self):
210 """ Comparison result, i.e. pass/warning/error """
211 self.ensure_compute()
212 return self._comparison_result
215 def comparison_result_long(self):
216 """ Longer description of the comparison result """
217 self.ensure_compute()
218 return self._comparison_result_long
222 """ This method performs the actual computations. """
225 def can_compare(self):
227 @return: True if the two objects can be compared, False otherwise
229 return self._has_correct_types()
and self._has_compatible_bins()
231 def _has_correct_types(self) -> bool:
233 @return: True if the two objects have a) a type supported for
234 comparison and b) can be compared with each other
236 if self.object_a
is None or self.object_b
is None:
241 supported_types = [
"TProfile",
"TH1D",
"TH1F",
"TEfficiency"]
242 if self.object_a.ClassName() != self.object_b.ClassName():
244 if self.object_a.ClassName()
not in supported_types:
247 if self.object_a.ClassName() ==
"TEfficiency":
249 if self.object_a.GetDimension() > 1:
254 def _raise_has_correct_types(self) -> None:
256 Raise Exception if not the two objects have a) a type supported for
257 comparison and b) can be compared with each other
260 if not self._has_correct_types():
261 msg =
"Comparison of {} (Type {}) with {} (Type {}) not " \
262 "supported.\nPlease open a JIRA issue (validation " \
263 "component) if you need this supported. "
264 raise ObjectsNotSupported(
266 self.object_a.GetName(),
267 self.object_a.ClassName(),
268 self.object_b.GetName(),
269 self.object_b.ClassName()
273 def _has_compatible_bins(self) -> bool:
275 Check if both ROOT obeject have the same amount of bins
276 @return: True if the bins are equal, otherwise False
278 if self.object_a.ClassName() ==
"TEfficiency" == \
279 self.object_b.ClassName():
280 nbins_a = self.object_a.GetTotalHistogram().GetNbinsX()
281 nbins_b = self.object_b.GetTotalHistogram().GetNbinsX()
283 nbins_a = self.object_a.GetNbinsX()
284 nbins_b = self.object_b.GetNbinsX()
286 return nbins_a == nbins_b
288 def _raise_has_compatible_bins(self) -> None:
290 Raise Exception if not both ROOT obeject have the same amount of bins
293 if not self._has_compatible_bins():
294 msg =
"The objects have differing x bin count: {} has {} vs. {} " \
296 raise DifferingBinCount(
298 self.object_a.GetName(),
299 self.object_a.GetNbinsX(),
300 self.object_b.GetName(),
301 self.object_b.GetNbinsX()
306 def _convert_teff_to_hist(teff_a):
308 Convert the content of a TEfficiency plot to a histogram and set
309 the bin content and errors
311 conv_hist = teff_a.GetTotalHistogram()
312 xbin_count = conv_hist.GetNbinsX()
313 xbin_low = conv_hist.GetXaxis().GetXmin()
314 xbin_max = conv_hist.GetXaxis().GetXmax()
317 teff_a.GetName() +
"root_conversion",
325 for i
in range(1, xbin_count):
326 th1.SetBinContent(i, teff_a.GetEfficiency(i))
327 th1.SetBinError(i, teff_a.GetEfficiencyErrorLow(i))
332 class PvalueTest(ComparisonBase):
333 """ Test with a pvalue """
337 _default_pvalue_warn = 1.0
341 _default_pvalue_error = 0.01
343 def __init__(self, *args, **kwargs):
344 """ Initialize Pvalue test
347 *args: Positional arguments to ComparisonBase
348 **kwargs: Keyword arguments to ComparisonBase
350 super().__init__(*args, **kwargs)
354 self._pvalue_warn = self.mop.pvalue_warn()
356 self._pvalue_error = self.mop.pvalue_error()
358 if self._pvalue_warn
is None:
359 self._pvalue_warn = self._default_pvalue_warn
360 if self._pvalue_error
is None:
361 self._pvalue_error = self._default_pvalue_error
363 def _get_comparison_result(self) -> str:
364 if self._pvalue
is None:
367 if self._pvalue < self._pvalue_error:
369 elif self._pvalue < self._pvalue_warn:
379 def _get_comparison_result_long(self):
391 class Chi2Test(PvalueTest):
394 Perform a Chi2Test for ROOT objects. The chi2 test method is e.g. described
395 in the documentation of TH1::Chi2Test. Basically this class wraps around
396 this Chi2Test function, and takes care that we can call perform these
397 tests for a wider selection of ROOT objects.
400 def __init__(self, *args, **kwargs):
403 :param args: See arguments of :class:`ComparisonBase`
404 :param kwargs: See arguments of :class:`ComparisonBase`
406 super().__init__(*args, **kwargs)
417 def _ensure_zero_error_has_no_content(self, a, b):
419 Ensure there are no bins which have a content set, but 0 error
420 This bin content will be set to 0 to disable this bin completely during
423 nbins = a.GetNbinsX()
424 for ibin
in range(1, nbins + 1):
425 if a.GetBinError(ibin) <= 0.0
and b.GetBinError(ibin) <= 0.0:
428 a.SetBinContent(ibin, 0.0)
429 b.SetBinContent(ibin, 0.0)
431 print(
"DEBUG: Warning: Setting bin content of bin {} to "
432 "zero for both histograms, because both histograms "
433 "have vanishing errors there.".format(ibin))
435 def _compute(self) -> None:
437 Performs the actual Chi^2 test
440 self._raise_has_correct_types()
441 self._raise_has_compatible_bins()
443 local_object_a = self.object_a
444 local_object_b = self.object_b
447 if self.object_a.ClassName() ==
"TEfficiency":
448 local_object_a = self._convert_teff_to_hist(self.object_a)
449 local_object_b = self._convert_teff_to_hist(self.object_b)
451 print(
"Converting TEfficiency objects to histograms.")
453 nbins = local_object_a.GetNbinsX()
456 raise TooFewBins(
"{} bin(s) is too few to perform the Chi2 "
457 "test.".format(nbins))
459 weighted_types = [
"TProfile",
"TH1D",
"TH1F"]
460 comp_weight_a = local_object_a.ClassName()
in weighted_types
461 comp_weight_b = local_object_b.ClassName()
in weighted_types
465 first_obj = local_object_a.Clone()
466 second_obj = local_object_b.Clone()
468 if comp_weight_a
and not comp_weight_b:
471 first_obj, second_obj = second_obj, first_obj
473 print(
"Debug: Warning: Switching the two objects, because "
474 "ROOT can only have the first one to be unweighted")
478 if comp_weight_a
and comp_weight_b:
480 elif comp_weight_a
or comp_weight_b:
485 if comp_weight_a
and comp_weight_b:
486 self._ensure_zero_error_has_no_content(first_obj, second_obj)
489 res_chi2 = numpy.array([1], numpy.float64)
490 res_igood = numpy.array([1], numpy.int32)
491 res_ndf = numpy.array([1], numpy.int32)
493 res_pvalue = first_obj.Chi2TestX(
502 print(
"Performing our own chi2 test, with bin-by-bin results: ")
504 print_contents_and_errors(first_obj, second_obj)
506 print(
"Here's what ROOT's Chi2Test gave us (comp_options: '{}'):"
507 " ".format(comp_options))
509 tp = TablePrinter(3, width=(10, 10, 40))
512 tp.print([
"Key",
"Value",
"Comment"])
514 tp.print([
"chi2", numpy.asscalar(res_chi2),
515 "Should roughly match above 'Total chi2'"])
516 tp.print([
"ndf", numpy.asscalar(res_ndf),
"#Non-empty bins - 1"])
517 tp.print([
"chi2/ndf", numpy.asscalar(res_chi2 / res_ndf),
""])
518 tp.print([
"igood", numpy.asscalar(res_igood),
519 "a debug indicator, 0 if all good"])
520 tp.print([
"pvalue", res_pvalue,
""])
523 print(
"See https://root.cern.ch/doc/master/classTH1.html for more "
528 msg =
"Comparison failed, no Chi^2 could be computed. For " \
529 "debugging, you can use the CLI of " \
530 "'validation/scripts/validationcomparison.py' on your root " \
531 "file and the reference. Run 'validationcomparison.py " \
532 "--help' for info. If problem persists, please open " \
533 "JIRA issue (validation component)."
534 raise ComparisonFailed(msg)
536 res_chi2ndf = res_chi2 / res_ndf
538 self._pvalue, self._chi2, self._chi2ndf, self._ndf = \
539 res_pvalue, res_chi2[0], res_chi2ndf[0], res_ndf[0]
541 def _get_comparison_result_long(self) -> str:
542 if self._pvalue
is None or self._chi2ndf
is None or self._chi2
is None:
543 return r"Could not perform $\chi^2$-Test between {{revision1}} " \
544 r"and {{revision2}} due to an unknown error. Please " \
545 r"submit a bug report."
547 return r'Performed $\chi^2$-Test between {{revision1}} ' \
548 r'and {{revision2}} ' \
549 r'($\chi^2$ = {chi2:.4f}; NDF = {ndf}; ' \
550 r'$\chi^2/\text{{{{NDF}}}}$ = {chi2ndf:.4f}).' \
551 r' <b>p-value: {pvalue:.6f}</b> (p-value warn: {pvalue_warn}, ' \
552 r'p-value error: {pvalue_error})'.format(
553 chi2=self._chi2, ndf=self._ndf, chi2ndf=self._chi2ndf,
554 pvalue=self._pvalue, pvalue_warn=self._pvalue_warn,
555 pvalue_error=self._pvalue_error
563 class KolmogorovTest(PvalueTest):
564 """ Kolmogorov-Smirnov Test """
566 def __init__(self, *args, **kwargs):
568 Initialize Kolmogorov test.
569 @param args: See arguments of :class:`ComparisonBase`
570 @param kwargs: See arguments of :class:`ComparisonBase`
572 super().__init__(*args, **kwargs)
576 Perform the actual test
579 self._raise_has_correct_types()
580 self._raise_has_compatible_bins()
582 local_object_a = self.object_a
583 local_object_b = self.object_b
586 if self.object_a.ClassName() ==
"TEfficiency":
587 local_object_a = self._convert_teff_to_hist(self.object_a)
588 local_object_b = self._convert_teff_to_hist(self.object_b)
590 print(
"Converting TEfficiency objects to histograms.")
596 self._pvalue = local_object_a.KolmogorovTest(local_object_b, option_str)
598 def _get_comparison_result_long(self) -> str:
599 if self._pvalue
is None:
600 return r"Could not perform Kolmogorov test between {{revision1}} " \
601 r"and {{revision2}} due to an unknown error. Please submit " \
604 return r'Performed Komlogorov test between {{revision1}} ' \
605 r'and {{revision2}} ' \
606 r' <b>p-value: {pvalue:.6f}</b> (p-value warn: {pvalue_warn}, ' \
607 r'p-value error: {pvalue_error})'.format(
608 pvalue=self._pvalue, pvalue_warn=self._pvalue_warn,
609 pvalue_error=self._pvalue_error
617 class AndersonDarlingTest(PvalueTest):
618 """ Anderson-Darling test"""
620 def __init__(self, *args, **kwargs):
622 Initialize Kolmogorov test.
623 @param args: See arguments of :class:`ComparisonBase`
624 @param kwargs: See arguments of :class:`ComparisonBase`
626 super().__init__(*args, **kwargs)
630 Perform the actual test
633 self._raise_has_correct_types()
638 local_object_a = self.object_a
639 local_object_b = self.object_b
642 if self.object_a.ClassName() ==
"TEfficiency":
643 local_object_a = self._convert_teff_to_hist(self.object_a)
644 local_object_b = self._convert_teff_to_hist(self.object_b)
646 print(
"Converting TEfficiency objects to histograms.")
652 self._pvalue = local_object_a.KolmogorovTest(local_object_b, option_str)
654 def _get_comparison_result_long(self) -> str:
655 if self._pvalue
is None:
656 return r"Could not perform Anderson Darling test between " \
657 r"{{revision1}} and {{revision2}} due to an unknown error." \
658 r" Please support a bug report."
660 return r'Performed Anderson Darling test between {{revision1}} ' \
661 r'and {{revision2}} ' \
662 r' <b>p-value: {pvalue:.6f}</b> (p-value warn: {pvalue_warn}, ' \
663 r'p-value error: {pvalue_error})'.format(
664 pvalue=self._pvalue, pvalue_warn=self._pvalue_warn,
665 pvalue_error=self._pvalue_error
673 class TablePrinter(object):
674 """ A tiny class to print columns of fixed width numbers. """
676 def __init__(self, ncols, width=None):
679 @param ncols: Number of columns
680 @param width: Width of each column. Either int or list.
686 if isinstance(width, int):
688 self.widths = [width] * ncols
689 elif isinstance(width, list)
or isinstance(width, tuple):
695 """ Total width of the table """
698 width += sum(self.widths)
700 width += (self.ncols - 1) * 3
705 def print_divider(self, char="="):
706 """ Print a divider made up from repeated chars """
707 print(char * self.tot_width)
709 def print(self, cols):
710 """ Print one row """
711 assert(len(cols) == self.ncols)
713 for icol, col
in enumerate(cols):
714 width = self.widths[icol]
715 if isinstance(col, int):
716 form =
"{{:{}d}}".format(width)
717 out.append(form.format(col))
718 elif isinstance(col, float):
719 form =
"{{:{}.{}f}}".format(width, width // 2)
720 out.append(form.format(col))
724 col = col[:width].rjust(width)
726 print(
"| " +
" | ".join(out) +
" |")
729 def print_contents_and_errors(obj_a, obj_b):
731 Print contents, errors and chi2 deviation for each bin as well as
732 some other information about two TH1-like objects.
733 @param obj_a: First TH1-like object
734 @param obj_b: Second TH1-like object
737 nbins = obj_a.GetNbinsX()
739 total_a = sum([obj_a.GetBinContent(ibin)
for ibin
in range(0, nbins + 2)])
740 total_b = sum([obj_b.GetBinContent(ibin)
for ibin
in range(0, nbins + 2)])
742 print(f
"Total events/summed weights in object 1: {total_a:10.5f}")
743 print(f
"Total events/summed weights in object 2: {total_b:10.5f}")
750 cp.print([
"ibin",
"a",
"err a",
"b",
"err b",
"chi2"])
752 for ibin
in range(1, nbins + 1):
753 content_a = obj_a.GetBinContent(ibin)
754 content_b = obj_b.GetBinContent(ibin)
755 error_a = obj_a.GetBinError(ibin)
756 error_b = obj_b.GetBinError(ibin)
760 chi2 = (total_b * content_a - total_a * content_b)**2 / \
761 (total_b**2 * error_a**2 + total_a**2 * error_b**2)
763 except ZeroDivisionError:
776 print(f
"Total chi2: {chi2_tot:10.5f}")
784 """ A small command line interface for debugging purposes. """
789 desc =
"For testing purposes: Run the chi2 comparison with objects from " \
791 parser = argparse.ArgumentParser(desc)
793 _ =
"Rootfile to read the first object from"
794 parser.add_argument(
"rootfile_a", help=_)
796 _ =
"Name of object inside first rootfile."
797 parser.add_argument(
"name_a", help=_)
799 _ =
"Rootfile to read the second object from"
800 parser.add_argument(
"rootfile_b", help=_)
802 _ =
"Name of object inside second rootfile."
803 parser.add_argument(
"name_b", help=_)
805 args = parser.parse_args()
810 if not os.path.exists(args.rootfile_a):
811 raise ValueError(f
"Could not find '{args.rootfile_a}'.")
813 if not os.path.exists(args.rootfile_b):
814 raise ValueError(f
"Could not find '{args.rootfile_b}'.")
816 rootfile_a = ROOT.TFile(args.rootfile_a)
817 obj_a = rootfile_a.Get(args.name_a)
820 f
"Could not find object '{args.name_a}' "
821 f
"in file '{args.rootfile_a}'.")
823 rootfile_b = ROOT.TFile(args.rootfile_b)
824 obj_b = rootfile_b.Get(args.name_b)
827 f
"Could not find object '{args.name_b}' "
828 f
"in file '{args.rootfile_b}'.")
833 test = Chi2Test(obj_a, obj_b, debug=
True)
834 test.ensure_compute()
836 print(
"If you see this message, then no exception was thrown.")
845 if __name__ ==
"__main__":