11""" Compare ROOT objects and perform e.g. chi2 tests.
12A small command line interface for testing/debugging purposes is included.
13Run `python3 validationcomparison.py --help` for more information. """
16from abc
import ABC, abstractmethod
20from typing
import Optional, Any
26from metaoptions
import MetaOptionParser
38class ComparisonFailed(Exception):
40 The comparison failed for some reason. For example
41 because ROOT was
not able to compute the Chi^2 properly
45class ObjectsNotSupported(Exception):
47 The type and/
or combination of provided ROOT objects
48 is not supported
for comparison
52class DifferingBinCount(Exception):
54 The two ROOT objects provided have a different bin count
55 and therefore, cannot be compared using the Chi2 test
59class TooFewBins(Exception):
61 Not sufficient bins to perform the Chi^2 test
65# ==============================================================================
66# Comparison class selector
67# ==============================================================================
71 object_1, object_2, mop: Optional[MetaOptionParser] = None,
73 """ Uses the metaoptions to determine which comparison algorithm is used
74 and initializes the corresponding subclass of :
class:`ComparisonBase` that
75 implements the actual comparison
and holds the results.
76 @param object_1 ROOT TObject
77 @param object_2 ROOT TObject
78 @param mop Metaoption parser
81 mop = MetaOptionParser()
82 if mop.has_option(
"kolmogorov"):
83 tester: Any = KolmogorovTest
84 elif mop.has_option(
"andersondarling"):
85 tester = AndersonDarlingTest
89 test = tester(object_1, object_2, mop=mop)
99class ComparisonBase(ABC):
101 Base class for all comparison implementations.
105 1. Initialize the
class together with two
ROOT objects of different revisions (that are to be compared)
and the metaoptions (given
in the
106 corresponding validation (steering) file), that determine how to compare
110 internally, but does
not compute anything yet
112 3. If :meth:`ensure_compute`
is called,
or any property
is accessed that
113 depends on computation, the internal implementation :meth:`_compute`
114 (to be implemented
in the subclass)
is called.
116 4. :meth:`_compute` ensures that all values, like chi2, p-value etc. are
119 5. Two properties :meth:`comparison_result` (
pass/warning/error)
and
120 :meth:`comparison_result_long` (longer description of the comparison result)
121 allow to access the results.
128 mop: Optional[MetaOptionParser] =
None,
132 Initialize ComparisonBase class
136 :param mop: MetaOptionParser
137 :param debug (bool): Debug mode enabled?
140 self.object_a = object_a
143 self.object_b = object_b
147 mop = MetaOptionParser()
154 self.computed =
False
157 self._comparison_result =
"not_compared"
160 self._comparison_result_long =
""
162 def ensure_compute(self):
164 Ensure all required quantities get computed and are cached inside the
170 if self.mop.has_option(
"nocompare"):
172 self._comparison_result_long =
"Testing is disabled for this plot"
175 fail_message =
"Comparison failed: "
180 except ObjectsNotSupported
as e:
181 self._comparison_result_long = fail_message + str(e)
182 except DifferingBinCount
as e:
183 self._comparison_result =
"error"
184 self._comparison_result_long = fail_message + str(e)
185 except TooFewBins
as e:
186 self._comparison_result_long = fail_message + str(e)
187 except ComparisonFailed
as e:
188 self._comparison_result =
"error"
189 self._comparison_result_long = fail_message + str(e)
190 except Exception
as e:
191 self._comparison_result =
"error"
192 self._comparison_result_long = (
193 "Unknown error occurred. Please "
194 "submit a bug report. " + str(e)
199 self._comparison_result_long = self._get_comparison_result_long()
200 self._comparison_result = self._get_comparison_result()
205 def _get_comparison_result(self) -> str:
206 """ Used to format the value of :attr:`_comparison_result`. """
209 def _get_comparison_result_long(self) -> str:
210 """ Used to format the value of :attr:`_comparison_result_long`. """
213 def comparison_result(self):
214 """ Comparison result, i.e. pass/warning/error """
215 self.ensure_compute()
216 return self._comparison_result
219 def comparison_result_long(self):
220 """ Longer description of the comparison result """
221 self.ensure_compute()
222 return self._comparison_result_long
226 """ This method performs the actual computations. """
228 def can_compare(self):
230 @return:
True if the two objects can be compared,
False otherwise
232 return self._has_correct_types()
and self._has_compatible_bins()
234 def _has_correct_types(self) -> bool:
236 @return:
True if the two objects have a) a type supported
for
237 comparison
and b) can be compared
with each other
239 if self.object_a
is None or self.object_b
is None:
244 supported_types = [
"TProfile",
"TH1D",
"TH1F",
"TEfficiency"]
245 if self.object_a.ClassName() != self.object_b.ClassName():
247 if self.object_a.ClassName()
not in supported_types:
250 if self.object_a.ClassName() ==
"TEfficiency":
252 if self.object_a.GetDimension() > 1:
257 def _raise_has_correct_types(self) -> None:
259 Raise Exception if not the two objects have a) a type supported
for
260 comparison
and b) can be compared
with each other
263 if not self._has_correct_types():
265 "Comparison of {} (Type {}) with {} (Type {}) not "
266 "supported.\nPlease open a GitLab issue (validation "
267 "label) if you need this supported. "
269 raise ObjectsNotSupported(
271 self.object_a.GetName(),
272 self.object_a.ClassName(),
273 self.object_b.GetName(),
274 self.object_b.ClassName(),
278 def _has_compatible_bins(self) -> bool:
280 Check if both ROOT objects have the same amount of bins
281 @return:
True if the bins are equal, otherwise
False
284 self.object_a.ClassName()
286 == self.object_b.ClassName()
288 nbins_a = self.object_a.GetTotalHistogram().GetNbinsX()
289 nbins_b = self.object_b.GetTotalHistogram().GetNbinsX()
291 nbins_a = self.object_a.GetNbinsX()
292 nbins_b = self.object_b.GetNbinsX()
294 return nbins_a == nbins_b
296 def _raise_has_compatible_bins(self) -> None:
298 Raise Exception if not both ROOT objects have the same amount of bins
301 if not self._has_compatible_bins():
303 "The objects have differing x bin count: {} has {} vs. {} "
306 raise DifferingBinCount(
308 self.object_a.GetName(),
309 self.object_a.GetNbinsX(),
310 self.object_b.GetName(),
311 self.object_b.GetNbinsX(),
316 def _convert_teff_to_hist(teff_a):
318 Convert the content of a TEfficiency plot to a histogram and set
319 the bin content
and errors
321 conv_hist = teff_a.GetTotalHistogram()
322 xbin_count = conv_hist.GetNbinsX()
323 xbin_low = conv_hist.GetXaxis().GetXmin()
324 xbin_max = conv_hist.GetXaxis().GetXmax()
327 teff_a.GetName() + "root_conversion",
335 for i
in range(1, xbin_count):
336 th1.SetBinContent(i, teff_a.GetEfficiency(i))
337 th1.SetBinError(i, teff_a.GetEfficiencyErrorLow(i))
342class PvalueTest(ComparisonBase):
343 """ Test with a pvalue """
347 _default_pvalue_warn = 0.99
351 _default_pvalue_error = 0.01
353 def __init__(self, *args, **kwargs):
354 """ Initialize Pvalue test
357 *args: Positional arguments to ComparisonBase
358 **kwargs: Keyword arguments to ComparisonBase
360 super().__init__(*args, **kwargs)
364 self._pvalue_warn = self.mop.pvalue_warn()
366 self._pvalue_error = self.mop.pvalue_error()
368 if self._pvalue_warn
is None:
369 self._pvalue_warn = self._default_pvalue_warn
370 if self._pvalue_error
is None:
371 self._pvalue_error = self._default_pvalue_error
373 def _get_comparison_result(self) -> str:
374 if self._pvalue
is None:
377 if self._pvalue < self._pvalue_error:
379 elif self._pvalue < self._pvalue_warn:
389 def _get_comparison_result_long(self):
402class Chi2Test(PvalueTest):
405 Perform a Chi2Test for ROOT objects. The chi2 test method
is e.g. described
406 in the documentation of TH1::Chi2Test. Basically this
class wraps around this Chi2Test function,
and takes care that we can call perform these
407 tests
for a wider selection of ROOT objects.
410 def __init__(self, *args, **kwargs):
413 :param args: See arguments of :class:`ComparisonBase`
414 :param kwargs: See arguments of :
class:`ComparisonBase`
416 super().__init__(*args, **kwargs)
427 def _ensure_zero_error_has_no_content(self, a, b):
429 Ensure there are no bins which have a content set, but 0 error
430 This bin content will be set to 0 to disable this bin completely during
433 nbins = a.GetNbinsX()
434 for ibin
in range(1, nbins + 1):
435 if a.GetBinError(ibin) <= 0.0
and b.GetBinError(ibin) <= 0.0:
438 a.SetBinContent(ibin, 0.0)
439 b.SetBinContent(ibin, 0.0)
442 f
"DEBUG: Warning: Setting bin content of bin {ibin} to zero for both histograms, because both " +
443 "histograms have vanishing errors there.")
445 def _compute(self) -> None:
447 Performs the actual Chi^2 test
450 self._raise_has_correct_types()
451 self._raise_has_compatible_bins()
453 local_object_a = self.object_a
454 local_object_b = self.object_b
457 if self.object_a.ClassName() ==
"TEfficiency":
458 local_object_a = self._convert_teff_to_hist(self.object_a)
459 local_object_b = self._convert_teff_to_hist(self.object_b)
461 print(
"Converting TEfficiency objects to histograms.")
463 nbins = local_object_a.GetNbinsX()
467 f
"{nbins} bin(s) is too few to perform the Chi2 test."
470 weighted_types = [
"TProfile",
"TH1D",
"TH1F"]
471 comp_weight_a = local_object_a.ClassName()
in weighted_types
472 comp_weight_b = local_object_b.ClassName()
in weighted_types
476 first_obj = local_object_a.Clone()
477 second_obj = local_object_b.Clone()
479 if comp_weight_a
and not comp_weight_b:
482 first_obj, second_obj = second_obj, first_obj
485 "Debug: Warning: Switching the two objects, because "
486 "ROOT can only have the first one to be unweighted"
491 if comp_weight_a
and comp_weight_b:
493 elif comp_weight_a
or comp_weight_b:
498 if comp_weight_a
and comp_weight_b:
499 self._ensure_zero_error_has_no_content(first_obj, second_obj)
502 res_chi2 = numpy.array([1], numpy.float64)
503 res_igood = numpy.array([1], numpy.int32)
504 res_ndf = numpy.array([1], numpy.int32)
506 res_pvalue = first_obj.Chi2TestX(
507 second_obj, res_chi2, res_ndf, res_igood, comp_options
511 print(
"Performing our own chi2 test, with bin-by-bin results: ")
513 print_contents_and_errors(first_obj, second_obj)
516 f
"Here's what ROOT's Chi2Test gave us (comp_options: '{comp_options}'): "
519 tp = TablePrinter(3, width=(10, 10, 40))
522 tp.print([
"Key",
"Value",
"Comment"])
527 numpy.asscalar(res_chi2),
528 "Should roughly match above 'Total chi2'",
531 tp.print([
"ndf", numpy.asscalar(res_ndf),
"#Non-empty bins - 1"])
532 tp.print([
"chi2/ndf", numpy.asscalar(res_chi2 / res_ndf),
""])
536 numpy.asscalar(res_igood),
537 "a debug indicator, 0 if all good",
540 tp.print([
"pvalue", res_pvalue,
""])
544 "See https://root.cern.ch/doc/master/classTH1.html for more "
551 "Comparison failed, no Chi^2 could be computed. For "
552 "debugging, you can use the CLI of "
553 "'validation/scripts/validationcomparison.py' on your root "
554 "file and the reference. Run 'validationcomparison.py "
555 "--help' for info. If problem persists, please open "
556 "GitLab issue (validation label)."
558 raise ComparisonFailed(msg)
560 res_chi2ndf = res_chi2 / res_ndf
562 self._pvalue, self._chi2, self._chi2ndf, self._ndf = (
569 def _get_comparison_result_long(self) -> str:
570 if self._pvalue
is None or self._chi2ndf
is None or self._chi2
is None:
572 r"Could not perform $\chi^2$-Test between {{revision1}} "
573 r"and {{revision2}} due to an unknown error. Please "
574 r"submit a bug report."
578 r"Performed $\chi^2$-Test between {{revision1}} "
579 r"and {{revision2}} "
580 r"($\chi^2$ = {chi2:.4f}; NDF = {ndf}; "
581 r"$\chi^2/\text{{{{NDF}}}}$ = {chi2ndf:.4f})."
582 r" <b>p-value: {pvalue:.6f}</b> (p-value warn: {pvalue_warn}, "
583 r"p-value error: {pvalue_error})".format(
586 chi2ndf=self._chi2ndf,
588 pvalue_warn=self._pvalue_warn,
589 pvalue_error=self._pvalue_error,
599class KolmogorovTest(PvalueTest):
600 """ Kolmogorov-Smirnov Test """
602 def __init__(self, *args, **kwargs):
604 Initialize Kolmogorov test.
605 @param args: See arguments of :
class:`ComparisonBase`
606 @param kwargs: See arguments of :
class:`ComparisonBase`
608 super().__init__(*args, **kwargs)
612 Perform the actual test
615 self._raise_has_correct_types()
616 self._raise_has_compatible_bins()
618 local_object_a = self.object_a
619 local_object_b = self.object_b
622 if self.object_a.ClassName() ==
"TEfficiency":
623 local_object_a = self._convert_teff_to_hist(self.object_a)
624 local_object_b = self._convert_teff_to_hist(self.object_b)
626 print(
"Converting TEfficiency objects to histograms.")
632 self._pvalue = local_object_a.KolmogorovTest(local_object_b, option_str)
634 def _get_comparison_result_long(self) -> str:
635 if self._pvalue
is None:
637 r"Could not perform Kolmogorov test between {{revision1}} "
638 r"and {{revision2}} due to an unknown error. Please submit "
643 r"Performed Kolmogorov test between {{revision1}} "
644 r"and {{revision2}} "
645 r" <b>p-value: {pvalue:.6f}</b> (p-value warn: {pvalue_warn}, "
646 r"p-value error: {pvalue_error})".format(
648 pvalue_warn=self._pvalue_warn,
649 pvalue_error=self._pvalue_error,
659class AndersonDarlingTest(PvalueTest):
660 """ Anderson-Darling test"""
662 def __init__(self, *args, **kwargs):
664 Initialize Anderson-Darling test.
665 @param args: See arguments of :
class:`ComparisonBase`
666 @param kwargs: See arguments of :
class:`ComparisonBase`
668 super().__init__(*args, **kwargs)
672 Perform the actual test
675 self._raise_has_correct_types()
680 local_object_a = self.object_a
681 local_object_b = self.object_b
684 if self.object_a.ClassName() ==
"TEfficiency":
685 local_object_a = self._convert_teff_to_hist(self.object_a)
686 local_object_b = self._convert_teff_to_hist(self.object_b)
688 print(
"Converting TEfficiency objects to histograms.")
694 self._pvalue = local_object_a.AndersonDarlingTest(
695 local_object_b, option_str
698 def _get_comparison_result_long(self) -> str:
699 if self._pvalue
is None:
701 r"Could not perform-Anderson Darling test between "
702 r"{{revision1}} and {{revision2}} due to an unknown error."
703 r" Please support a bug report."
707 r"Performed Anderson-Darling test between {{revision1}} "
708 r"and {{revision2}} "
709 r" <b>p-value: {pvalue:.6f}</b> (p-value warn: {pvalue_warn}, "
710 r"p-value error: {pvalue_error})".format(
712 pvalue_warn=self._pvalue_warn,
713 pvalue_error=self._pvalue_error,
724 """ A tiny class to print columns of fixed width numbers. """
726 def __init__(self, ncols, width=None):
729 @param ncols: Number of columns
730 @param width: Width of each column. Either int
or list.
736 if isinstance(width, int):
738 self.widths = [width] * ncols
739 elif isinstance(width, list)
or isinstance(width, tuple):
745 """ Total width of the table """
748 width += sum(self.widths)
750 width += (self.ncols - 1) * 3
755 def print_divider(self, char="="):
756 """ Print a divider made up from repeated chars """
757 print(char * self.tot_width)
759 def print(self, cols):
760 """ Print one row """
761 assert len(cols) == self.ncols
763 for icol, col
in enumerate(cols):
764 width = self.widths[icol]
765 if isinstance(col, int):
766 form = f
"{{:{width}d}}"
767 out.append(form.format(col))
768 elif isinstance(col, float):
769 form = f
"{{:{width}.{width // 2}f}}"
770 out.append(form.format(col))
774 col = col[:width].rjust(width)
776 print(
"| " +
" | ".join(out) +
" |")
779def print_contents_and_errors(obj_a, obj_b):
781 Print contents, errors and chi2 deviation
for each bin
as well
as
782 some other information about two TH1-like objects.
783 @param obj_a: First TH1-like object
784 @param obj_b: Second TH1-like object
787 nbins = obj_a.GetNbinsX()
789 total_a = sum([obj_a.GetBinContent(ibin) for ibin
in range(0, nbins + 2)])
790 total_b = sum([obj_b.GetBinContent(ibin)
for ibin
in range(0, nbins + 2)])
792 print(f
"Total events/summed weights in object 1: {total_a:10.5f}")
793 print(f
"Total events/summed weights in object 2: {total_b:10.5f}")
800 cp.print([
"ibin",
"a",
"err a",
"b",
"err b",
"chi2"])
802 for ibin
in range(1, nbins + 1):
803 content_a = obj_a.GetBinContent(ibin)
804 content_b = obj_b.GetBinContent(ibin)
805 error_a = obj_a.GetBinError(ibin)
806 error_b = obj_b.GetBinError(ibin)
810 chi2 = (total_b * content_a - total_a * content_b) ** 2 / (
811 total_b ** 2 * error_a ** 2 + total_a ** 2 * error_b ** 2
814 except ZeroDivisionError:
816 cp.print([ibin, content_a, error_a, content_b, error_b, chi2])
820 print(f
"Total chi2: {chi2_tot:10.5f}")
829 """ A small command line interface for debugging purposes. """
835 "For testing purposes: Run the chi2 comparison with objects from "
838 parser = argparse.ArgumentParser(desc)
840 _ =
"Rootfile to read the first object from"
841 parser.add_argument(
"rootfile_a", help=_)
843 _ =
"Name of object inside first rootfile."
844 parser.add_argument(
"name_a", help=_)
846 _ =
"Rootfile to read the second object from"
847 parser.add_argument(
"rootfile_b", help=_)
849 _ =
"Name of object inside second rootfile."
850 parser.add_argument(
"name_b", help=_)
852 args = parser.parse_args()
857 if not os.path.exists(args.rootfile_a):
858 raise ValueError(f
"Could not find '{args.rootfile_a}'.")
860 if not os.path.exists(args.rootfile_b):
861 raise ValueError(f
"Could not find '{args.rootfile_b}'.")
863 rootfile_a = ROOT.TFile(args.rootfile_a)
864 obj_a = rootfile_a.Get(args.name_a)
867 f
"Could not find object '{args.name_a}' "
868 f
"in file '{args.rootfile_a}'."
871 rootfile_b = ROOT.TFile(args.rootfile_b)
872 obj_b = rootfile_b.Get(args.name_b)
875 f
"Could not find object '{args.name_b}' "
876 f
"in file '{args.rootfile_b}'."
882 test = Chi2Test(obj_a, obj_b, debug=
True)
883 test.ensure_compute()
885 print(
"If you see this message, then no exception was thrown.")
894if __name__ ==
"__main__":