17 from .
import statistics
23 """Get for the logging.Logger instance of this module
28 Logger instance of this module
30 return logging.getLogger(__name__)
34 units_by_quantity_name = {
59 def get_unit(quantity_name):
60 """Infers the unit of a quantity from its name.
62 Assumes the standard Belle II unit system.
64 Currently looks up the quantity string from units_by_quantity_name.
69 Name of a quantity (E.g. pt, x, ...)
76 unit = units_by_quantity_name.get(quantity_name,
None)
80 def compose_axis_label(quantity_name, unit=None):
81 """Formats a quantity name and a unit to a label for a plot axes.
83 If the unit is not given to is tried to infer it
84 from the quantity name by the get_unit function.
89 Name of the quantity to be displayed at the axes
91 The unit of the quantity. Defaults to get_unit(quantity_name)
99 unit = get_unit(quantity_name)
102 axis_label = quantity_name
104 axis_label =
'%s (%s)' % (quantity_name, unit)
109 def get1DBinningFromReference(name, refFileName):
110 """ returns nbins, lowerbound, upperbound for TH1 / TProfile with name "name" found in the file "refFileName"
112 @param name : name of the TH1 object to be looked for in the file
113 @param refFileName : name of the reference file where the object is searched for
115 @return int nbin, float xmin, float xmax of the TH1
122 if refFileName
is None or refFileName ==
"":
123 return nbins, x_min, x_max
126 oldDirectory = ROOT.gROOT.CurrentDirectory()
128 tfile = ROOT.TFile(refFileName)
130 objptr = tfile.Get(name)
131 if objptr
and objptr.InheritsFrom(
"TH1"):
132 nbins = objptr.GetNbinsX()
133 x_min = objptr.GetXaxis().GetXmin()
134 x_max = objptr.GetXaxis().GetXmax()
136 basf2.B2WARNING(
'Requested object with name: ' + name +
' not found in file: ' + refFileName +
" (or not a TH1)")
138 basf2.B2WARNING(
'Requested file: ' + refFileName +
' could not be opened')
145 return nbins, x_min, x_max
149 StatsEntry = ROOT.TParameter(float)
154 """Class for generating a validation plot for the Belle II validation page.
156 Typically it generates plots from values stored in numpy arrays and feeds them into
157 plot ROOT classes for storing them.
159 It implements an automatic binning procedure based on the rice rule and
160 robust z score outlier detection.
162 It also keeps track of additional statistics typically neglected by ROOT such as a count
163 for the non finit values such as NaN, +Inf, -Inf.
165 The special attributes for the Belle II validation page like
170 are exposed as properties of this class.
174 very_sparse_dots_line_style_index = 28
178 """Constructor of the ValidationPlot
183 A unique name to be used as the name of the ROOT object to be generated
185 referenceFileName : str
186 name of a reference file. If set the code will try to get the histogram or profile
187 from that file and determine the number of bins and upper and lower bound
188 (so far only implemented for 1D (TH1, TProfile), is ignored for 2D plots)
192 self.
name = root_save_name(name)
240 outlier_z_score=None,
241 include_exceptionals=True,
242 allow_discrete=False,
243 cumulation_direction=None,
245 """Fill the plot with a one dimensional histogram."""
250 if n
is not None and xmin
is not None and xmax
is not None:
255 th1_factory = ROOT.TH1D
263 lower_bound=lower_bound,
264 upper_bound=upper_bound,
265 outlier_z_score=outlier_z_score,
266 include_exceptionals=include_exceptionals,
267 allow_discrete=allow_discrete,
268 cumulation_direction=cumulation_direction)
286 outlier_z_score=None,
287 include_exceptionals=True,
288 allow_discrete=False,
289 cumulation_direction=None,
292 """Fill the plot with a one dimensional profile of one variable over another."""
300 if n
is not None and xmin
is not None and xmax
is not None:
305 th1_factory = ROOT.TProfile
307 if gaus_z_score
is None:
314 lower_bound=lower_bound,
315 upper_bound=upper_bound,
316 outlier_z_score=outlier_z_score,
317 include_exceptionals=include_exceptionals,
318 allow_discrete=allow_discrete,
319 cumulation_direction=cumulation_direction)
323 self.
hist2d(xs, ys=ys, weights=weights, stackby=stackby,
325 lower_bound=(lower_bound,
None),
326 upper_bound=(upper_bound,
None),
327 outlier_z_score=(outlier_z_score, outlier_z_score),
328 include_exceptionals=(include_exceptionals,
True),
329 allow_discrete=(allow_discrete,
False),
336 name=histogram.GetName()[1:],
337 z_score=gaus_z_score)
338 profiles.append(profile)
347 self.
ylabel =
'probability'
350 histogram.SetMinimum(0)
351 histogram.SetMaximum(1.05)
353 self.
plot.SetMinimum(0)
354 self.
plot.SetMaximum(1.05)
362 lower_bound=(
None,
None),
363 upper_bound=(
None,
None),
364 outlier_z_score=(
None,
None),
365 include_exceptionals=(
True,
True),
368 """Fill the plot with a (unbinned) two dimensional scatter plot"""
374 x_outlier_z_score, y_outlier_z_score = self.
unpack_2d_param(outlier_z_score)
375 x_include_exceptionals, y_include_exceptionals = self.
unpack_2d_param(include_exceptionals)
379 lower_bound=x_lower_bound,
380 upper_bound=x_upper_bound,
381 outlier_z_score=x_outlier_z_score,
382 include_exceptionals=x_include_exceptionals
387 lower_bound=y_lower_bound,
388 upper_bound=y_upper_bound,
389 outlier_z_score=y_outlier_z_score,
390 include_exceptionals=y_include_exceptionals
393 graph = ROOT.TGraph()
395 graph.SetName(self.
name)
396 graph.SetMarkerStyle(6)
397 graph.GetHistogram().SetOption(
"AP")
402 graph.SetLineColorAlpha(color_index, 0)
406 graph.GetXaxis().SetLimits(x_lower_bound, x_upper_bound)
407 graph.GetYaxis().SetLimits(y_lower_bound, y_upper_bound)
421 lower_bound=(
None,
None),
422 upper_bound=(
None,
None),
423 outlier_z_score=(
None,
None),
424 include_exceptionals=(
True,
True),
427 """Fill the plot with a (unbinned) two dimensional scatter plot
428 xs_and_err and ys_and_err are tuples containing the values and the errors on these values
439 x_outlier_z_score, y_outlier_z_score = self.
unpack_2d_param(outlier_z_score)
440 x_include_exceptionals, y_include_exceptionals = self.
unpack_2d_param(include_exceptionals)
444 lower_bound=x_lower_bound,
445 upper_bound=x_upper_bound,
446 outlier_z_score=x_outlier_z_score,
447 include_exceptionals=x_include_exceptionals
452 lower_bound=y_lower_bound,
453 upper_bound=y_upper_bound,
454 outlier_z_score=y_outlier_z_score,
455 include_exceptionals=y_include_exceptionals
458 graph = ROOT.TGraphErrors()
460 graph.SetName(self.
name)
461 graph.GetHistogram().SetOption(
"A")
463 graph.SetMarkerColor(4)
464 graph.SetMarkerStyle(21)
467 graph.GetXaxis().SetLimits(x_lower_bound, x_upper_bound)
468 graph.GetYaxis().SetLimits(y_lower_bound, y_upper_bound)
484 lower_bound=(
None,
None),
485 upper_bound=(
None,
None),
486 outlier_z_score=(
None,
None),
487 include_exceptionals=(
True,
True),
488 allow_discrete=(
False,
False),
491 """Fill the plot with a two dimensional histogram"""
495 if quantiles
is not None:
496 name =
"_" + self.
name
503 x_outlier_z_score, y_outlier_z_score = self.
unpack_2d_param(outlier_z_score)
504 x_include_exceptionals, y_include_exceptionals = self.
unpack_2d_param(include_exceptionals)
505 x_allow_discrete, y_allow_discrete = self.
unpack_2d_param(allow_discrete)
507 if quantiles
is not None:
508 y_include_exceptionals =
True
509 y_allow_discrete =
False
514 lower_bound=x_lower_bound,
515 upper_bound=x_upper_bound,
516 outlier_z_score=x_outlier_z_score,
517 include_exceptionals=x_include_exceptionals,
518 allow_discrete=x_allow_discrete)
523 lower_bound=y_lower_bound,
524 upper_bound=y_upper_bound,
525 outlier_z_score=y_outlier_z_score,
526 include_exceptionals=y_include_exceptionals,
527 allow_discrete=y_allow_discrete)
529 n_x_bins = len(x_bin_edges) - 1
530 n_y_bins = len(y_bin_edges) - 1
532 self.
lower_bound = (x_bin_edges[0], y_bin_edges[0])
534 self.
upper_bound = (x_bin_edges[-1], y_bin_edges[-1])
536 histogram = ROOT.TH2D(name,
544 get_logger().info(
"Scatter plot %s is discrete in x.", name)
545 x_taxis = histogram.GetXaxis()
546 for i_x_bin, x_bin_label
in enumerate(x_bin_labels):
547 x_taxis.SetBinLabel(i_x_bin + 1, x_bin_label)
551 x_bin_width = x_bin_edges[1] - x_bin_edges[0]
555 get_logger().info(
"Scatter plot %s is discrete in y.", name)
556 y_taxis = histogram.GetYaxis()
557 for i_y_bin, y_bin_label
in enumerate(y_bin_labels):
558 y_taxis.SetBinLabel(i_y_bin + 1, y_bin_label)
562 y_bin_width = y_bin_edges[1] - y_bin_edges[0]
565 self.
create(histogram, xs, ys=ys, weights=weights, stackby=stackby)
567 if quantiles
is not None:
571 for quantile
in quantiles:
572 profile = histogram.QuantilesX(quantile, histogram.GetName()[1:] +
'_' + str(quantile))
575 x_taxis = histogram.GetXaxis()
576 new_x_taxis = profile.GetXaxis()
577 for i_bin
in range(x_taxis.GetNbins() + 2):
578 label = x_taxis.GetBinLabel(i_bin)
580 new_x_taxis.SetBinLabel(i_bin, label)
583 epsilon = sys.float_info.epsilon
584 for i_bin
in range(0, profile.GetNbinsX() + 2):
585 profile.SetBinError(i_bin, epsilon)
587 profiles.append(profile)
590 self.
plot = self.
create_stack(profiles, name=self.
plot.GetName()[1:], reverse_stack=
False, force_graph=
True)
595 x_taxis = histogram.GetXaxis()
596 x_bin_edges = array.array(
"d", list(range(len(x_bin_labels) + 1)))
597 x_taxis.Set(n_x_bins, x_bin_edges)
602 x_taxis = histogram.GetXaxis()
603 y_bin_edges = array.array(
"d", list(range(len(y_bin_labels) + 1)))
604 y_taxis.Set(n_y_bins, y_bin_edges)
609 """Fit a gaus belle curve to the central portion of a one dimensional histogram
611 The fit is applied to the central mean +- z_score * std interval of the histogram,
612 such that it is less influence by non gaussian tails further away than the given z score.
614 @param float z_score number of sigmas to include from the mean value of the histogram.
621 raise RuntimeError(
'Validation plot must be filled before it can be fitted.')
623 if not isinstance(plot, ROOT.TH1D):
624 raise RuntimeError(
'Fitting is currently implemented / tested for one dimensional, non stacked validation plots.')
628 fit_tf1 = ROOT.TF1(
"Fit", formula)
629 fit_tf1.SetTitle(title)
630 fit_tf1.SetParName(0,
"n")
631 fit_tf1.SetParName(1,
"mean")
632 fit_tf1.SetParName(2,
"std")
634 n = histogram.GetSumOfWeights()
635 mean = histogram.GetMean()
636 std = histogram.GetStdDev()
638 fit_tf1.SetParameter(0, n)
639 fit_tf1.SetParameter(1, mean)
640 fit_tf1.SetParameter(2, std)
643 return self.
fit(fit_tf1,
648 """Fit a general line to a one dimensional histogram"""
651 fit_tf1 = ROOT.TF1(
"Fit", formula)
652 fit_tf1.SetTitle(title)
653 fit_tf1.SetParName(0,
"slope")
654 fit_tf1.SetParName(1,
"intercept")
655 self.
fit(fit_tf1,
'M')
658 """Fit a constant function to a one dimensional histogram"""
661 fit_tf1 = ROOT.TF1(
"Fit", formula)
662 fit_tf1.SetTitle(title)
663 fit_tf1.SetParName(0,
"intercept")
664 self.
fit(fit_tf1,
'M')
667 """Fit a diagonal line through the origin to a one dimensional histogram"""
670 fit_tf1 = ROOT.TF1(
"Fit", formula)
671 fit_tf1.SetTitle(title)
672 fit_tf1.SetParName(0,
"slope")
673 self.
fit(fit_tf1,
'M')
675 def fit(self, formula, options, lower_bound=None, upper_bound=None, z_score=None):
676 """Fit a user defined function to a one dimensional histogram
681 Formula string or TH1 to be fitted. See TF1 constructurs for that is a valid formula
683 Options string to be used in the fit. See TH1::Fit()
685 Lower bound of the range to be fitted
687 Upper bound of the range to be fitted
691 raise RuntimeError(
'Validation plot must be filled before it can be fitted.')
693 if not isinstance(plot, ROOT.TH1D):
694 raise RuntimeError(
'Fitting is currently implemented / tested for one dimensional, non stacked validation plots.')
698 xaxis = histogram.GetXaxis()
699 n_bins = xaxis.GetNbins()
700 hist_lower_bound = xaxis.GetBinLowEdge(1)
701 hist_upper_bound = xaxis.GetBinUpEdge(n_bins)
703 if z_score
is not None:
704 mean = histogram.GetMean()
705 std = histogram.GetStdDev()
707 if lower_bound
is None:
708 lower_bound = mean - z_score * std
710 if upper_bound
is None:
711 upper_bound = mean + z_score * std
714 if isinstance(formula, ROOT.TF1):
716 fit_tf1.SetRange(hist_lower_bound, hist_upper_bound)
718 fit_tf1 = ROOT.TF1(
"Fit",
722 get_logger().info(
'Fitting with %s', fit_tf1.GetExpFormula())
725 if lower_bound
is None or lower_bound < hist_lower_bound:
726 lower_bound = hist_lower_bound
727 if upper_bound
is None or upper_bound > hist_upper_bound:
728 upper_bound = hist_upper_bound
732 if 'N' not in options:
735 fit_res = histogram.Fit(fit_tf1, options +
"S",
"", lower_bound, upper_bound)
745 raise ValueError(
"Can not show a validation plot that has not been filled.")
748 """Write the plot to file
752 tdirectory : ROOT.TDirectory, optional
753 ROOT directory to which the plot should be written.
754 If omitted write to the current directory
757 raise ValueError(
"Can not write a validation plot that has not been filled.")
759 with root_cd(tdirectory):
760 ValidationPlot.set_tstyle()
765 meta_options = [
"nostats"]
769 meta_options.append(
"expert")
771 meta_options.append(
"shifter")
775 meta_options.append(
"pvalue-error={}".format(self.
pvalue_error))
777 meta_options.append(
"pvalue-warn={}".format(self.
pvalue_warn))
781 meta_options.append(
"logy")
783 meta_options_str =
",".join(meta_options)
786 histogram.GetListOfFunctions().Add(ROOT.TNamed(
'MetaOptions', meta_options_str))
791 """Getter method if an plot plot is marked as expert plot"""
796 """Getter for the plot title"""
801 """Setter for the plot title"""
804 self.
plot.SetTitle(title)
806 histogram.SetTitle(title)
810 """Getter for the axis label at the x axis"""
815 """Setter for the axis label at the x axis"""
818 histogram.GetXaxis().SetTitle(xlabel)
822 """Getter for the axis label at the y axis"""
827 """Setter for the axis label at the y axis"""
830 histogram.GetYaxis().SetTitle(ylabel)
834 """Getter for the contact email address to be displayed on the validation page"""
839 """Setter for the contact email address to be displayed on the validation page"""
842 found_obj = histogram.FindObject(
'Contact')
844 tnamed = ROOT.TNamed(
"Contact", contact)
845 histogram.GetListOfFunctions().Add(tnamed)
846 found_obj = histogram.FindObject(
'Contact')
847 found_obj.SetTitle(contact)
851 """Getter for the description to be displayed on the validation page"""
856 """Setter for the description to be displayed on the validation page"""
859 found_obj = histogram.FindObject(
'Description')
861 tnamed = ROOT.TNamed(
"Description", description)
862 histogram.GetListOfFunctions().Add(tnamed)
863 found_obj = histogram.FindObject(
'Description')
864 found_obj.SetTitle(description)
868 """Getter for the check to be displayed on the validation page"""
873 """Setter for the check to be displayed on the validation page"""
876 found_obj = histogram.FindObject(
'Check')
878 tnamed = ROOT.TNamed(
"Check", check)
879 histogram.GetListOfFunctions().Add(tnamed)
880 found_obj = histogram.FindObject(
'Check')
881 found_obj.SetTitle(check)
888 """Unpacks a function parameter for the two dimensional plots.
890 If it is a pair the first parameter shall apply to the x coordinate
891 the second to the y coordinate. In this case the pair is returned as two values
893 If something else is given the it is assumed that this parameter should equally apply
894 to both coordinates. In this case the same values is return twice as a pair.
898 param : pair or single value
899 Function parameter for a two dimensional plot
904 A pair of values being the parameter for the x coordinate and
905 the y coordinate respectively
909 x_param, y_param = param
913 return x_param, y_param
917 """Determine if the data consists of boolean values"""
918 return statistics.is_binary_series(xs)
922 """Determine if the data consists of discrete values"""
923 return statistics.is_discrete_series(xs, max_n_unique=max_n_unique)
927 """Find exceptionally frequent values
937 A list of the found exceptional values.
939 return statistics.rice_exceptional_values(xs)
943 """Does an estimation of mean and standard deviation robust against outliers.
953 Pair of mean and standard deviation
955 x_mean = statistics.truncated_mean(xs)
956 x_std = statistics.trimmed_std(xs)
961 """Formats a value to be placed at a tick on an axis."""
962 if np.isfinite(value)
and value == np.round(value):
963 return str(int(value))
965 formated_value =
"{:.5g}".format(value)
968 if len(formated_value) > 8:
969 formated_value =
"{:.3e}".format(value)
970 return formated_value
981 outlier_z_score=None,
982 include_exceptionals=True,
983 allow_discrete=False,
984 cumulation_direction=None):
985 """Combined factory method for creating a one dimensional histogram or a profile plot."""
989 xs = np.array(xs, copy=
False)
992 ys = np.array(ys, copy=
False)
994 if weights
is not None:
995 weights = np.array(weights, copy=
False)
1000 lower_bound=lower_bound,
1001 upper_bound=upper_bound,
1002 outlier_z_score=outlier_z_score,
1003 include_exceptionals=include_exceptionals,
1004 allow_discrete=allow_discrete)
1006 n_bins = len(bin_edges) - 1
1009 histogram = th1_factory(name,
'', n_bins, bin_edges)
1012 get_logger().info(
"One dimensional plot %s is discrete in x.", name)
1013 x_taxis = histogram.GetXaxis()
1014 for i_bin, bin_label
in enumerate(bin_labels):
1015 x_taxis.SetBinLabel(i_bin + 1, bin_label)
1019 bin_width = bin_edges[1] - bin_edges[0]
1027 cumulation_direction=cumulation_direction,
1034 x_taxis = histogram.GetXaxis()
1035 bin_edges = array.array(
"d", list(range(len(bin_labels) + 1)))
1036 x_taxis.Set(n_bins, bin_edges)
1044 cumulation_direction=None,
1045 reverse_stack=None):
1046 """Create histograms from a template, possibly stacked"""
1051 histogram = histogram_template
1052 self.
fill_into(histogram, xs, ys, weights=weights)
1053 if cumulation_direction
is not None:
1054 histogram = self.
cumulate(histogram, cumulation_direction=cumulation_direction)
1056 histograms.append(histogram)
1060 stackby = np.array(stackby, copy=
False)
1061 name = histogram_template.GetName()
1068 groupby_label=
"stack")
1070 if cumulation_direction
is not None:
1071 histograms = [self.
cumulate(histogram, cumulation_direction=cumulation_direction)
1072 for histogram
in histograms]
1074 plot = self.
create_stack(histograms, name=name +
"_stacked", reverse_stack=reverse_stack)
1082 """Create a stack of histograms"""
1083 if len(histograms) == 1:
1084 plot = histograms[0]
1086 if isinstance(histograms[0], (ROOT.TProfile, ROOT.TGraph))
or force_graph:
1087 plot = ROOT.TMultiGraph()
1089 plot = ROOT.THStack()
1096 for histogram
in reversed(histograms):
1097 if isinstance(histogram, ROOT.TProfile)
or (isinstance(histogram, ROOT.TH1)
and force_graph):
1099 plot.Add(histogram,
"APZ")
1103 for histogram
in histograms:
1104 if isinstance(histogram, ROOT.TProfile)
or (isinstance(histogram, ROOT.TH1)
and force_graph):
1106 plot.Add(histogram,
"APZ")
1114 """Extract errors from a TProfile histogram and create a TGraph from these"""
1115 if isinstance(tprofile, ROOT.TGraph):
1118 x_taxis = tprofile.GetXaxis()
1119 n_bins = x_taxis.GetNbins()
1121 bin_ids_with_underflow = list(range(n_bins + 1))
1122 bin_ids_without_underflow = list(range(1, n_bins + 1))
1124 bin_centers = np.array([x_taxis.GetBinCenter(i_bin)
for i_bin
in bin_ids_without_underflow])
1126 bin_centers = np.abs(bin_centers)
1127 bin_widths = np.array([x_taxis.GetBinWidth(i_bin)
for i_bin
in bin_ids_without_underflow])
1128 bin_x_errors = bin_widths / 2.0
1131 bin_contents = np.array([tprofile.GetBinContent(i_bin)
for i_bin
in bin_ids_without_underflow])
1132 bin_y_errors = np.array([tprofile.GetBinError(i_bin)
for i_bin
in bin_ids_without_underflow])
1134 tgrapherrors = ROOT.TGraphErrors(n_bins, bin_centers, bin_contents, bin_x_errors, bin_y_errors)
1136 tgrapherrors.GetHistogram().SetOption(
"APZ")
1138 tgrapherrors.SetLineColor(tprofile.GetLineColor())
1139 tgrapherrors.SetLineColor(tprofile.GetLineColor())
1142 for tobject
in tprofile.GetListOfFunctions():
1143 tgrapherrors.GetListOfFunctions().Add(tobject.Clone())
1148 stats_values = np.array([np.nan] * 6)
1149 tprofile.GetStats(stats_values)
1151 sum_w = stats_values[0]
1152 sum_w2 = stats_values[1]
1153 sum_wx = stats_values[2]
1154 sum_wx2 = stats_values[3]
1155 sum_wy = stats_values[4]
1156 sum_wy2 = stats_values[5]
1164 np.sqrt(sum_wx2 * sum_w - sum_wx * sum_wx) / sum_w)
1172 np.sqrt(sum_wy2 * sum_w - sum_wy * sum_wy) / sum_w)
1176 tgrapherrors.GetCovariance())
1180 tgrapherrors.GetCorrelationFactor())
1190 groupby_label="group"):
1191 """Fill data into similar histograms in groups indicated by a groupby array"""
1194 unique_groupbys = np.unique(groupbys)
1195 name = histogram_template.GetName()
1197 for i_value, value
in enumerate(unique_groupbys):
1199 indices_for_value = np.isnan(groupbys)
1201 indices_for_value = groupbys == value
1204 histogram_for_value = histogram_template.Clone(name +
'_' + str(value))
1205 i_root_color = i_value + 1
1207 self.
set_color(histogram_for_value, i_root_color)
1216 filter=indices_for_value)
1218 histograms.append(histogram_for_value)
1223 """Set the color of the ROOT object.
1225 By default the line color of a TGraph should be invisible, so do not change it
1226 For other objects set the marker and the line color
1230 tobject : Plotable object inheriting from TAttLine and TAttMarker such as TGraph or TH1
1231 Object of which the color should be set.
1233 Color index of the ROOT color table
1235 if isinstance(tobject, ROOT.TGraph):
1236 tobject.SetMarkerColor(root_i_color)
1238 tobject.SetLineColor(root_i_color)
1239 tobject.SetMarkerColor(root_i_color)
1241 def fill_into(self, plot, xs, ys=None, weights=None, filter=None):
1242 """Fill the data into the plot object"""
1243 if isinstance(plot, ROOT.TGraph):
1245 raise ValueError(
"ys are required for filling a graph")
1247 elif isinstance(plot, ROOT.TGraphErrors):
1249 raise ValueError(
"ys are required for filling a graph error")
1253 self.
fill_into_th1(plot, xs, ys, weights=weights, filter=filter)
1256 """fill point values and error of the x and y axis into the graph"""
1258 assert(len(xs[0]) == len(ys[0]))
1260 graph.Set(len(xs[0]))
1262 for i
in range(len(xs[0])):
1263 graph.SetPoint(i, xs[0][i], ys[0][i])
1264 graph.SetPointError(i, xs[1][i], ys[1][i])
1267 """Fill the data into a TGraph"""
1271 filter =
slice(
None)
1281 if x_n_data > max_n_data
or y_n_data > max_n_data:
1282 get_logger().warning(
"Number of points in scatter graph %s exceed limit %s" %
1283 (self.
name, max_n_data))
1285 get_logger().warning(
"Cropping %s" % max_n_data)
1287 xs = xs[0:max_n_data]
1288 ys = ys[0:max_n_data]
1290 x_axis = graph.GetXaxis()
1291 y_axis = graph.GetYaxis()
1293 x_lower_bound = x_axis.GetXmin()
1294 x_upper_bound = x_axis.GetXmax()
1296 y_lower_bound = y_axis.GetXmin()
1297 y_upper_bound = y_axis.GetXmax()
1299 x_underflow_indices = xs < x_lower_bound
1300 x_overflow_indices = xs > x_upper_bound
1302 y_underflow_indices = ys < y_lower_bound
1303 y_overflow_indices = ys > y_upper_bound
1305 plot_indices = ~(np.isnan(xs) |
1306 x_underflow_indices |
1307 x_overflow_indices |
1309 y_underflow_indices |
1312 n_plot_data = np.sum(plot_indices)
1313 plot_xs = xs[plot_indices]
1314 plot_ys = ys[plot_indices]
1316 graph.Set(int(n_plot_data))
1317 for i, (x, y)
in enumerate(zip(plot_xs, plot_ys)):
1318 graph.SetPoint(i, x, y)
1325 x_n_underflow = np.sum(x_underflow_indices)
1329 x_n_overflow = np.sum(x_overflow_indices)
1333 y_n_underflow = np.sum(y_underflow_indices)
1337 y_n_overflow = np.sum(y_overflow_indices)
1351 """Fill the histogram blocking non finite values
1355 histogram : ROOT.TH1
1356 The histogram to be filled
1357 xs : numpy.ndarray (1d)
1358 Data for the first axes
1359 ys : numpy.ndarray (1d), optional
1360 Data for the second axes
1361 weights : numpy.ndarray (1d), optional
1362 Weight of the individual points. Defaults to one for each
1363 filter : numpy.ndarray, optional
1364 Boolean index array indicating which entries shall be taken.
1368 filter =
slice(
None)
1373 finite_filter = np.isfinite(xs)
1379 finite_filter &= np.isfinite(ys)
1382 xs = xs[finite_filter]
1383 weights = np.ones_like(xs)
1385 weights = weights[filter]
1387 finite_filter &= np.isfinite(weights)
1388 xs = xs[finite_filter]
1389 weights[finite_filter]
1392 ys = ys[finite_filter]
1397 except AttributeError:
1398 Fill = histogram.Fill
1400 fill = np.frompyfunc(Fill, 2, 1)
1401 fill(xs.astype(np.float64, copy=
False),
1402 weights.astype(np.float64, copy=
False))
1404 fill = np.frompyfunc(Fill, 3, 1)
1405 fill(xs.astype(np.float64, copy=
False),
1406 ys.astype(np.float64, copy=
False),
1407 weights.astype(np.float64, copy=
False))
1411 xs = xs.astype(np.float64, copy=
False)
1412 weights = weights.astype(np.float64, copy=
False)
1415 histogram.FillN(n, xs, weights)
1417 basf2.B2WARNING(
"No values to be filled into histogram: " + self.
name)
1421 xs = xs.astype(np.float64, copy=
False)
1422 ys = ys.astype(np.float64, copy=
False)
1423 weights = weights.astype(np.float64, copy=
False)
1426 histogram.FillN(n, xs, ys, weights)
1428 basf2.B2WARNING(
"No values to be filled into histogram: " + self.
name)
1434 """ Extracts the counts of non finite floats from a series
1435 and adds them as additional statistics to the histogram.
1439 histogram : derived from ROOT.TH1 or ROOT.TGraph
1440 Something having a GetListOfFunctions method that
1442 A label for the data series to be prefixed to the entries.
1443 xs : numpy.ndarray (1d)
1444 Data from which the non finit floats should be counted.
1446 n_nans = np.isnan(xs).sum()
1450 n_positive_inf = np.sum(xs == np.inf)
1451 if n_positive_inf > 0:
1454 n_negative_inf = np.sum(xs == -np.inf)
1455 if n_negative_inf > 0:
1460 """Add a new additional statistics to the histogram.
1464 histogram : derived from ROOT.TH1 or ROOT.TGraph
1465 Something having a GetListOfFunctions method that holds the additional statistics
1467 Label of the statistic
1469 Value of the statistic
1471 stats_entry = StatsEntry(str(label), float(value))
1472 histogram.GetListOfFunctions().Add(stats_entry)
1477 """Get the additional statistics from the histogram and return them a dict.
1481 histogram : derived from ROOT.TH1 or ROOT.TGraph
1482 Something having a GetListOfFunctions method that holds the additional statistics
1486 collection.OrderedDict
1487 A map of labels to values for the additional statistics
1489 additional_stats = collections.OrderedDict()
1490 for tobject
in histogram.GetListOfFunctions():
1491 if isinstance(tobject, StatsEntry):
1492 stats_entry = tobject
1493 label = stats_entry.GetName()
1494 value = stats_entry.GetVal()
1495 additional_stats[label] = value
1496 return additional_stats
1500 """Extract a slice of a scatterplot and apply a Gaussian fit to it"""
1503 y_taxis = th2.GetYaxis()
1504 th2_lower_bound = y_taxis.GetXmin()
1505 th2_upper_bound = y_taxis.GetXmax()
1506 th2_height = y_taxis.GetXmax() - y_taxis.GetXmin()
1507 n_y_bins = y_taxis.GetNbins()
1509 y_mean = th2.GetMean(2)
1510 y_std = th2.GetStdDev(2)
1511 fit_lower_bound = max(th2_lower_bound, y_mean - z_score * y_std)
1512 fit_upper_bound = min(th2_upper_bound, y_mean + z_score * y_std)
1513 fit_height = fit_upper_bound - fit_lower_bound
1515 required_n_bins_inslice_filled = n_y_bins * fit_height / th2_height
1517 fit_lower_bound = th2_lower_bound
1518 fit_upper_bound = th2_upper_bound
1519 fit_height = fit_upper_bound - fit_lower_bound
1520 required_n_bins_inslice_filled = n_y_bins / 1.61
1523 required_n_bins_inslice_filled = min(required_n_bins_inslice_filled, n_y_bins / 1.61)
1525 fit_tf1 = ROOT.TF1(
"Fit",
"gaus", fit_lower_bound, fit_upper_bound)
1526 fit_tf1.SetParName(0,
"n")
1527 fit_tf1.SetParName(1,
"mean")
1528 fit_tf1.SetParName(2,
"std")
1532 param_fit_th1s = ROOT.TObjArray()
1533 th2.FitSlicesY(fit_tf1, i_first_bin, i_last_bin,
1534 int(required_n_bins_inslice_filled),
1535 fit_options, param_fit_th1s)
1537 th1_means = param_fit_th1s.At(1)
1538 th1_means.SetName(name)
1539 th1_means.SetTitle(th2.GetTitle())
1544 x_taxis = th2.GetXaxis()
1545 new_x_taxis = th1_means.GetXaxis()
1546 for i_bin
in range(x_taxis.GetNbins() + 2):
1547 label = x_taxis.GetBinLabel(i_bin)
1549 new_x_taxis.SetBinLabel(i_bin, label)
1552 data_lower_bound = th1_means.GetMinimum(fit_lower_bound)
1553 data_upper_bound = th1_means.GetMaximum(fit_upper_bound)
1554 data_height = data_upper_bound - data_lower_bound
1556 plot_lower_bound = max(fit_lower_bound, data_lower_bound - 0.05 * data_height)
1557 plot_upper_bound = min(fit_upper_bound, data_upper_bound + 0.05 * data_height)
1559 th1_means.SetMinimum(plot_lower_bound)
1560 th1_means.SetMaximum(plot_upper_bound)
1565 def cumulate(cls, histogram, cumulation_direction=None):
1566 """Cumulates the histogram inplace.
1570 histogram : ROOT.TH1 or ROOT.TProfile
1571 Filled histogram to be cumulated
1572 cumulation_direction : int, optional
1573 Direction is indicated by the sign.
1574 Positive means from left to right, negative means from right to left.
1575 If now cumulation direction is given return the histogram as is.
1580 Cumulated histogram potentially altered inplace.
1582 if not cumulation_direction:
1585 cumulate_backward = cumulation_direction < 0
1586 cumulate_forward =
not cumulate_backward
1588 if isinstance(histogram, ROOT.TH2):
1589 raise ValueError(
"Cannot cumulate a two dimensional histogram.")
1591 if isinstance(histogram, ROOT.TH3):
1592 raise ValueError(
"Cannot cumulate a three dimensional histogram.")
1594 if not isinstance(histogram, ROOT.TH1):
1595 raise ValueError(
"Can only cumulate a one dimensional histogram.")
1597 if isinstance(histogram, ROOT.TProfile):
1598 tprofile = histogram
1601 tgraph.SetName(tprofile.GetName())
1603 n_bins = histogram.GetNbinsX()
1605 cumulated_content = 0.0
1606 cumulated_entries = 0
1610 i_bins = list(range(0, n_bins + 2))
1611 if not cumulate_forward:
1612 i_bins = reversed(i_bins)
1614 for i_bin
in i_bins:
1616 bin_content = tprofile.GetBinContent(i_bin)
1617 bin_entries = tprofile.GetBinEffectiveEntries(i_bin)
1618 bin_std = tprofile.GetBinError(i_bin)
1620 if bin_entries != 0:
1621 cumulated_content = (
1622 1.0 * (cumulated_entries * cumulated_content + bin_entries * bin_content) /
1623 (cumulated_entries + bin_entries)
1627 math.hypot(cumulated_entries * cumulated_std, bin_entries * bin_std) /
1628 (cumulated_entries + bin_entries)
1631 cumulated_entries = cumulated_entries + bin_entries
1636 if i_point >= 0
and i_point < n_points:
1637 x = tgraph.GetX()[i_point]
1641 tgraph.SetPoint(i_point, x, cumulated_content)
1643 x_error = tgraph.GetErrorX(i_point)
1644 tgraph.SetPointError(i_point, x_error, cumulated_std)
1649 n_bins = histogram.GetNbinsX()
1650 cumulated_content = 0.0
1652 i_bins = list(range(0, n_bins + 2))
1653 if not cumulate_forward:
1654 i_bins = reversed(i_bins)
1656 for i_bin
in i_bins:
1657 bin_content = histogram.GetBinContent(i_bin)
1658 cumulated_content += bin_content
1659 histogram.SetBinContent(i_bin, cumulated_content)
1669 outlier_z_score=None,
1670 include_exceptionals=True,
1671 allow_discrete=False):
1672 """Deducing bin edges from a data series.
1676 xs : numpy.ndarray (1d)
1677 Data point for which a binning should be found.
1678 stackbys : numpy.ndarray (1d)
1679 Categories of the data points to be distinguishable
1680 bins : list(float) or int or None, optional
1681 Preset bin edges or preset number of desired bins.
1682 The default, None, means the bound should be extracted from data.
1683 The rice rule is used the determine the number of bins.
1684 If a list of floats is given return them immediatly.
1685 lower_bound : float or None, optional
1686 Preset lower bound of the binning range.
1687 The default, None, means the bound should be extracted from data.
1688 upper_bound : float or None, optional
1689 Preset upper bound of the binning range.
1690 The default, None, means the bound should be extracted from data.
1691 outlier_z_score : float or None, optional
1692 Threshold z-score of outlier detection.
1693 The default, None, means no outlier detection.
1694 include_exceptionals : bool, optional
1695 If the outlier detection is active this switch indicates,
1696 if values detected as exceptionally frequent shall be included
1697 nevertheless into the binning range. Default is True,
1698 which means exceptionally frequent values as included
1699 even if they are detected as outliers.
1703 np.array (1d), list(str)
1704 Pair of bin edges and labels deduced from the series.
1705 Second element is None if the series is not detected as discrete.
1707 debug = get_logger().debug
1708 debug(
'Determine binning for plot named %s', self.
name)
1714 elif isinstance(bins, collections.Iterable):
1718 bin_edges = array.array(
'd', bin_edges)
1720 return bin_edges, bin_labels
1731 message =
'Cannot accept n_bins=%s as number of bins, because it is not a number greater than 0.' % bins
1732 raise ValueError(message)
1735 xs = np.array(xs, copy=
False)
1739 debug(
'Discrete binning values encountered')
1740 finite_xs = xs[np.isfinite(xs)]
1741 unique_xs = np.unique(finite_xs)
1744 if lower_bound
is None:
1745 if len(unique_xs) == 0:
1746 if upper_bound
is None:
1749 lower_bound = upper_bound - 1
1751 lower_bound = np.min(unique_xs)
1753 unique_xs = unique_xs[unique_xs >= lower_bound]
1755 if upper_bound
is None:
1756 if len(unique_xs) == 0:
1757 upper_bound = lower_bound + 1
1759 upper_bound = np.min(unique_xs)
1761 unique_xs = unique_xs[unique_xs <= upper_bound]
1764 n_bins = len(unique_xs)
or 1
1766 if len(unique_xs) > 0
and n_bins >= len(unique_xs):
1768 bin_edges = array.array(
'd', unique_xs)
1771 bin_edges.append(bin_edges[-1] + 1)
1772 return bin_edges, bin_labels
1781 debug(
'Lower bound %s', lower_bound)
1782 debug(
'Upper bound %s', upper_bound)
1783 debug(
'N bins %s', n_bins)
1789 lower_bound=lower_bound,
1790 upper_bound=upper_bound,
1791 outlier_z_score=outlier_z_score,
1792 include_exceptionals=include_exceptionals)
1794 n_bins, lower_bound, upper_bound = bin_range
1796 n_bin_edges = n_bins + 1
1797 if lower_bound != upper_bound:
1799 debug(
"Creating flat distribution binning")
1800 precentiles = np.linspace(0.0, 100.0, n_bin_edges)
1801 bin_edges = np.unique(np.nanpercentile(xs[(lower_bound <= xs) & (xs <= upper_bound)], precentiles))
1805 bin_edges = np.linspace(lower_bound, upper_bound, n_bin_edges)
1810 bin_edges[0] = lower_bound
1811 bin_edges[-1] = np.nextafter(upper_bound, np.inf)
1812 debug(
'Bins %s', bin_edges)
1816 bin_edges = [lower_bound, upper_bound + 1]
1819 bin_edges = array.array(
'd', bin_edges)
1820 debug(
'Bins %s for %s', bin_edges, self.
name)
1821 return bin_edges,
None
1829 outlier_z_score=None,
1830 include_exceptionals=True):
1831 """Calculates the number of bins, the lower bound and the upper bound from a given data series
1832 estimating the values that are not given.
1834 If the outlier_z_score is given the method tries to exclude outliers that exceed a certain z-score.
1835 The z-score is calculated (x - x_mean) / x_std. The be robust against outliers the necessary
1836 mean and std deviation are based on truncated mean and a trimmed std calculated from the inter
1837 quantile range (IQR).
1839 If additional include_exceptionals is true the method tries to find exceptional values in the series
1840 and always include them in the range if it finds any.
1841 Exceptional values means exact values that appear often in the series for whatever reason.
1842 Possible reasons include
1843 * Interal / default values
1844 * Failed evaluation conditions
1846 which should be not cropped away automatically if you are locking on the quality of your data.
1850 xs : numpy.ndarray (1d)
1851 Data point for which a binning should be found.
1852 stackbys : numpy.ndarray (1d)
1853 Categories of the data points to be distinguishable
1854 n_bins : int or None, optional
1855 Preset number of desired bins. The default, None, means the bound should be extracted from data.
1856 The rice rule is used the determine the number of bins.
1857 lower_bound : float or None, optional
1858 Preset lower bound of the binning range.
1859 The default, None, means the bound should be extracted from data.
1860 upper_bound : float or None, optional
1861 Preset upper bound of the binning range.
1862 The default, None, means the bound should be extracted from data.
1863 outlier_z_score : float or None, optional
1864 Threshold z-score of outlier detection.
1865 The default, None, means no outlier detection.
1866 include_exceptionals : bool, optional
1867 If the outlier detection is active this switch indicates,
1868 if values detected as exceptionally frequent shall be included
1869 nevertheless into the binning range. Default is True,
1870 which means exceptionally frequent values as included
1871 even if they are detected as outliers.
1875 n_bins, lower_bound, upper_bound : int, float, float
1876 A triple of found number of bins, lower bound and upper bound of the binning range.
1879 if stackbys
is not None:
1880 unique_stackbys = np.unique(stackbys)
1882 for value
in unique_stackbys:
1884 indices_for_value = np.isnan(stackbys)
1886 indices_for_value = stackbys == value
1888 stack_lower_bound, stack_upper_bound = \
1890 lower_bound=lower_bound,
1891 upper_bound=upper_bound,
1892 outlier_z_score=outlier_z_score,
1893 include_exceptionals=include_exceptionals)
1895 stack_ranges.append([stack_lower_bound, stack_upper_bound])
1897 lower_bound = np.nanmin([lwb
for lwb, upb
in stack_ranges])
1898 upper_bound = np.nanmax([upb
for lwb, upb
in stack_ranges])
1902 lower_bound=lower_bound,
1903 upper_bound=upper_bound,
1904 outlier_z_score=outlier_z_score,
1905 include_exceptionals=include_exceptionals)
1910 n_data = np.sum((lower_bound <= xs) & (xs <= upper_bound))
1911 rice_n_bins = int(statistics.rice_n_bin(n_data))
1912 n_bins = rice_n_bins
1915 n_bins = int(n_bins)
1918 message =
'Cannot accept n_bins=%s as number of bins, because it is not a number greater than 0.' % bins
1919 raise ValueError(message)
1921 return n_bins, lower_bound, upper_bound
1927 outlier_z_score=None,
1928 include_exceptionals=True):
1932 xs : numpy.ndarray (1d)
1933 Data point for which a binning should be found.
1934 lower_bound : float or None, optional
1935 Preset lower bound of the binning range.
1936 The default, None, means the bound should be extracted from data.
1937 upper_bound : float or None, optional
1938 Preset upper bound of the binning range.
1939 The default, None, means the bound should be extracted from data.
1940 outlier_z_score : float or None, optional
1941 Threshold z-score of outlier detection.
1942 The default, None, means no outlier detection.
1943 include_exceptionals : bool, optional
1944 If the outlier detection is active this switch indicates,
1945 if values detected as exceptionally frequent shall be included
1946 nevertheless into the binning range. Default is True,
1947 which means exceptionally frequent values as included
1948 even if they are detected as outliers.
1952 lower_bound, upper_bound : float, float
1953 A pair of found lower bound and upper bound of series.
1955 debug = get_logger().debug
1957 finite_xs_indices = np.isfinite(xs)
1958 if np.any(finite_xs_indices):
1959 finite_xs = xs[finite_xs_indices]
1963 make_symmetric =
False
1964 exclude_outliers = outlier_z_score
is not None and (lower_bound
is None or upper_bound
is None)
1967 if include_exceptionals
or exclude_outliers:
1969 exceptional_indices = np.in1d(finite_xs, exceptional_xs)
1972 if exclude_outliers:
1973 if not np.all(exceptional_indices):
1978 x_mean, x_std = np.nan, np.nan
1980 make_symmetric = abs(x_mean) < x_std / 5.0
and lower_bound
is None and upper_bound
is None
1982 if include_exceptionals
and len(exceptional_xs) != 0:
1983 lower_exceptional_x = np.min(exceptional_xs)
1984 upper_exceptional_x = np.max(exceptional_xs)
1985 make_symmetric =
False
1987 lower_exceptional_x = np.nan
1988 upper_exceptional_x = np.nan
1991 if lower_bound
is None:
1993 lower_bound = np.min(finite_xs)
1997 if outlier_z_score
is not None:
1999 lower_outlier_bound = x_mean - outlier_z_score * x_std
2003 indices_above_lower_outlier_bound = finite_xs >= lower_outlier_bound
2005 if np.any(indices_above_lower_outlier_bound):
2006 lower_bound = np.min(finite_xs[indices_above_lower_outlier_bound])
2009 lower_bound = np.nanmin([lower_bound, lower_exceptional_x])
2011 debug(
'Lower bound after outlier detection')
2012 debug(
'Lower bound %s', lower_bound)
2013 debug(
'Lower outlier bound %s', lower_outlier_bound)
2016 if upper_bound
is None:
2018 upper_bound = np.max(finite_xs)
2021 if outlier_z_score
is not None:
2023 upper_outlier_bound = x_mean + outlier_z_score * x_std
2027 indices_below_upper_outlier_bound = finite_xs <= upper_outlier_bound
2029 if np.any(indices_below_upper_outlier_bound):
2030 upper_bound = np.max(finite_xs[indices_below_upper_outlier_bound])
2033 upper_bound = np.nanmax([upper_bound, upper_exceptional_x])
2035 debug(
'Upper bound after outlier detection')
2036 debug(
'Upper bound %s', upper_bound)
2037 debug(
'Upper outlier bound %s', upper_outlier_bound)
2039 if make_symmetric
and lower_bound < 0
and upper_bound > 0:
2040 if abs(abs(lower_bound) - abs(upper_bound)) < x_std / 5.0:
2041 abs_bound = max(abs(lower_bound), abs(upper_bound))
2042 lower_bound = -abs_bound
2043 upper_bound = abs_bound
2045 return lower_bound, upper_bound
2049 """Combining fit TF1 with the additional statistics and attach them to the histogram.
2053 histogram : ROOT.TH1 or ROOT.TGraph or ROOT.TMultiGraph
2054 Something having a GetListOfFunctions method that should hold
2055 the combined fit and additional statistics function.
2058 cls.
set_tf1(histogram, additional_stats_tf1)
2062 """Combining fit TF1 with the additional statistics and attach them to the histogram.
2066 histogram : ROOT.TH1 or ROOT.TGraph or ROOT.TMultiGraph
2067 Something having a GetListOfFunctions method that should hold
2068 the combined fit and additional statistics function.
2072 cls.
set_tf1(histogram, combined_tf1)
2076 """Set the attached TF1 of the histogram.
2080 histogram : ROOT.TH1 or ROOT.TGraph or ROOT.TMultiGraph
2081 Something having a GetListOfFunctions method that should hold
2082 the combined fit and additional statistics function.
2087 tf1.SetName(
"FitAndStats")
2088 histogram.GetListOfFunctions().Add(tf1)
2092 """Delete the attached TF1 from the histogram
2096 histogram : ROOT.TH1 or ROOT.TGraph
2097 Something having a GetListOfFunctions method that holds the fit function
2099 tf1 = histogram.FindObject(
"FitAndStats")
2101 function_list = histogram.GetListOfFunctions()
2102 function_list.Remove(tf1)
2106 """Create a TF1 with the additional statistics from the histogram as parameters.
2110 histogram : ROOT.TH1 or ROOT.TGraph
2111 Something having a GetListOfFunctions method that holds the additional statistics.
2116 Function with the additional statistics as parameters.
2120 if not additional_stats:
2130 formula_string =
'+'.join(
'0*[' + str(i) +
']' for i
in range(len(additional_stats)))
2133 additional_stats_tf1 = ROOT.TF1(
"Stats", formula_string, lower_bound, upper_bound)
2135 for (i, (label, value))
in enumerate(additional_stats.items()):
2139 label = label.replace(
" ",
"-")
2140 additional_stats_tf1.SetParName(i, label)
2141 additional_stats_tf1.FixParameter(i, value)
2143 return additional_stats_tf1
2147 """Combine the fit function and the function carrying the additional statistics to one function.
2153 additional_stats_tf1 : ROOT.TF1
2154 The function carrying the additional statistics as parameters
2160 if additional_stats_tf1
is None:
2166 lower_bound = ROOT.Double()
2167 upper_bound = ROOT.Double()
2168 fit_tf1.GetRange(lower_bound, upper_bound)
2169 title = fit_tf1.GetTitle()
2171 combined_formula = additional_stats_tf1.GetExpFormula().Data() +
'+' + fit_tf1.GetExpFormula().Data()
2172 combined_tf1 = ROOT.TF1(
"Combined", combined_formula, lower_bound, upper_bound)
2173 combined_tf1.SetTitle(title)
2176 chi2 = fit_tf1.GetChisquare()
2177 combined_tf1.SetChisquare(chi2)
2179 ndf = fit_tf1.GetNDF()
2180 combined_tf1.SetNDF(ndf)
2182 n_stats_parameters = additional_stats_tf1.GetNpar()
2183 n_fit_parameters = fit_tf1.GetNpar()
2191 """Copy the parameters of one TF1 to another.
2195 tf1_source : ROOT.TF1
2196 Take parameters from here
2197 tf1_target : ROOT.TF1
2199 offset : int, optional
2200 Index of the first target parameter to which to copy.
2202 n_parameters = tf1_source.GetNpar()
2205 lower_bound = ROOT.Double()
2206 upper_bound = ROOT.Double()
2208 for i_source
in range(n_parameters):
2209 parameter_name = tf1_source.GetParName(i_source)
2210 i_target = tf1_target.GetParNumber(parameter_name)
2214 for i_target
in range(target_tf1.GetNpar()):
2215 if parameter_name == target_tf1.GetParName(i_target):
2221 tf1_target.SetParameter(i_target,
2222 tf1_source.GetParameter(i_source))
2223 tf1_target.SetParError(i_target,
2224 tf1_source.GetParError(i_source))
2226 tf1_source.GetParLimits(i_source, lower_bound, upper_bound)
2227 tf1_target.SetParLimits(i_target, lower_bound, upper_bound)
2230 """Reassign the special attributes of the plot forwarding them to the ROOT plot."""
2247 """Sets the maximum of the vertical plotable range"""
2249 if isinstance(histogram, ROOT.TH1):
2250 histogram.SetMaximum(histogram.GetMaximum(maximum))
2252 histogram.SetMaximum(maximum)
2255 """Sets the minimum of the vertical plotable range"""
2257 if isinstance(histogram, ROOT.TH1):
2258 histogram.SetMinimum(histogram.GetMinimum(minimum))
2260 histogram.SetMinimum(minimum)
2264 """Set the style such that the additional stats entries are shown by the TBrowser"""
2265 belle2_validation_style_name =
"belle2_validation_style"
2266 belle2_validation_tstyle = ROOT.gROOT.GetStyle(belle2_validation_style_name)
2267 if not belle2_validation_tstyle:
2268 belle2_validation_tstyle = ROOT.TStyle(belle2_validation_style_name, belle2_validation_style_name)
2271 belle2_validation_tstyle.SetOptFit(opt_fit)
2274 belle2_validation_tstyle.SetOptStat(opt_stat)
2275 ROOT.gROOT.SetStyle(belle2_validation_style_name)
2280 belle2_validation_tstyle.cd()
2284 """Simple test methode"""
2285 ValidationPlot.set_tstyle()
2288 normal_distributed_values = np.random.randn(1000)
2291 normal_distributed_values[i] = np.nan
2293 for i
in range(10, 20):
2294 normal_distributed_values[i] = np.inf
2296 for i
in range(20, 30):
2297 normal_distributed_values[i] = -np.inf
2300 validation_histogram.hist(normal_distributed_values)
2301 validation_histogram.title =
'A normal distribution'
2302 validation_histogram.xlabel =
'normal'
2303 validation_histogram.ylabel =
'frequency'
2304 validation_histogram.fit_gaus()
2308 cumulated_histogram.hist(normal_distributed_values, cumulation_direction=1)
2309 cumulated_histogram.title =
'A cumulated normal distribution'
2310 cumulated_histogram.xlabel =
'normal'
2311 cumulated_histogram.ylabel =
'cdf'
2312 cumulated_histogram.show()
2316 stackby = np.random.binomial(1.0, 0.40, 1000)
2317 stacked_validation_histogram =
ValidationPlot(
'test_stacked_hist')
2318 stacked_validation_histogram.hist(normal_distributed_values, stackby=stackby)
2321 x = np.random.randn(1000)
2322 y = 3 * np.random.randn(1000)
2323 ones = np.ones_like(x)
2326 x1 = np.where(stackby != 0, np.cos(angle) * ones, ones) * x + np.where(stackby != 0, np.sin(angle) * ones, ones) * y
2327 y1 = np.where(stackby != 0, np.sin(angle) * ones, ones) * x - np.where(stackby != 0, np.cos(angle) * ones, ones) * y
2329 stacked_validation_scatter =
ValidationPlot(
'test_stacked_scatter')
2330 stacked_validation_scatter.scatter(x1, y1, stackby=stackby)
2333 stacked_validation_profile =
ValidationPlot(
'test_stacked_profile')
2334 stacked_validation_profile.profile(x1, y1, stackby=stackby)
2337 stacked_validation_hist2d =
ValidationPlot(
'test_stacked_hist2d')
2338 stacked_validation_hist2d.hist2d(x1, y1, stackby=stackby)
2341 x = np.linspace(-1, 1, 1000)
2345 diagonal_plot.profile(x, y, bins=50)
2346 diagonal_plot.fit_line()
2350 cumulated_profile.profile(x, y, bins=50, cumulation_direction=1)
2352 tfile = ROOT.TFile(
'test.root',
'RECREATE')
2354 validation_histogram.write(tfile)
2356 with root_cd(
"expert")
as tdirectory1:
2357 diagonal_plot.write(tdirectory1)
2358 cumulated_profile.write(tdirectory1)
2359 cumulated_histogram.write(tdirectory1)
2361 with root_cd(
"stacked")
as tdirectory2:
2362 stacked_validation_histogram.write(tdirectory2)
2363 stacked_validation_scatter.write()
2364 stacked_validation_profile.write()
2365 stacked_validation_hist2d.write()
2369 tfile = ROOT.TFile(
'test.root')
2370 tBrowser = ROOT.TBrowser()
2371 tBrowser.BrowseObject(tfile)
2376 if __name__ ==
'__main__':