Belle II Software  release-05-01-25
metadata.py
1 #!/usr/bin/env python3
2 
3 """Small module containing helper functions to set the metadata on objects
4 created for the validation correctly """
5 
6 # std
7 from typing import Optional, Union, List, Tuple
8 import pathlib
9 
10 import basf2
11 import ROOT
12 from ROOT import Belle2
13 
14 # circumvent BII-1264
15 ROOT.gInterpreter.Declare("#include <framework/utilities/MakeROOTCompatible.h>")
16 
17 # Unfortunately doxygen doesn't recognize our docstrings in this file :/
18 # @cond SUPPRESS_DOXYGEN
19 
20 
21 def file_description_set(rootfile: Union[ROOT.TFile, str, pathlib.PurePath],
22  description: str) -> None:
23  """
24  Add file description validation metdata to a ROOT file.
25 
26  Args:
27  rootfile (TFile, str or pathlib.PurePath): Name of the root file
28  to open or an already open TFile instance
29  description (str): Common description/information of/about all plots
30  in this ROOT file (will be displayed above the plots)
31 
32  Returns:
33  None
34  """
35  opened = False
36  if not isinstance(rootfile, ROOT.TFile):
37  if isinstance(rootfile, pathlib.PurePath):
38  rootfile = str(rootfile)
39  rootfile = ROOT.TFile(rootfile, "UPDATE")
40  opened = True
41  if not rootfile.IsOpen() or not rootfile.IsWritable():
42  raise RuntimeError(
43  f"ROOT file {rootfile.GetName()} is not open for writing"
44  )
45  # scope guard to avoid side effects by changing the global gDirectory
46  # in modules ...
47  # noinspection PyUnusedLocal
48  directory_guard = ROOT.TDirectory.TContext(rootfile)
49  desc = ROOT.TNamed("Description", description)
50  desc.Write()
51  if opened:
52  rootfile.Close()
53 
54 
55 def validation_metadata_set(obj: ROOT.TObject, title: str, contact: str,
56  description: str, check: str,
57  xlabel: Optional[str] = None,
58  ylabel: Optional[str] = None,
59  metaoptions="") -> None:
60  """
61  Set the validation metadata for a given object by setting the necessary
62  values. This function can be used on any object supported by the
63  Validation (histograms, profiles, ntuples)
64 
65  Arguments:
66  obj: Instance of the object which should get the metadata
67  title (str): Title to use for the object
68  contact (str): Contact person, usually in the form "Name <email>"
69  description (str): Text description what can be seen in the plot
70  check (str): Text description what to look for in the validation for
71  shifters to easily see if the distribution looks ok
72  xlabel (str): If given try to set this as the label for the x axis
73  ylabel (str): If given try to set this as the label for the y axis
74  metaoptions (str): Metaoptions (additional options to influence the
75  comparison between revisions, styling of the plot, etc.)
76 
77  .. warning::
78 
79  Different ways to specify LaTeX for different arguments:
80  see `create_validation_histograms`
81  """
82  obj.SetTitle(title)
83  # For ntuples we add the metadata as aliases ...
84  try:
85  obj.SetAlias("Contact", contact)
86  obj.SetAlias("Description", description)
87  obj.SetAlias("Check", check)
88  obj.SetAlias("MetaOptions", metaoptions)
89  except AttributeError:
90  pass
91  # for TH*, TProfile we add it to the list of functions
92  try:
93  function_list = obj.GetListOfFunctions()
94  function_list.Add(ROOT.TNamed("Contact", contact))
95  function_list.Add(ROOT.TNamed("Description", description))
96  function_list.Add(ROOT.TNamed("Check", check))
97  function_list.Add(ROOT.TNamed("MetaOptions", metaoptions))
98  except AttributeError:
99  pass
100 
101  # And maybe try to set a label on the axis
102  if xlabel is not None:
103  try:
104  obj.GetXaxis().SetTitle(xlabel)
105  except AttributeError:
106  pass
107 
108  if ylabel is not None:
109  try:
110  obj.GetYaxis().SetTitle(ylabel)
111  except AttributeError:
112  pass
113 
114 
115 # noinspection PyIncorrectDocstring
116 def validation_metadata_update(
117  rootfile: Union[str, ROOT.TFile, pathlib.PurePath],
118  name: str, *args, **argk) -> None:
119  """
120  This is a convenience helper for `validation_metadata_set` in case the
121  objects have already been saved in a ROOT file before: It will open the
122  file (or use an existing TFile), extract the object, add the metadata and
123  save the new version to the file
124 
125  Arguments:
126  rootfile (str, ROOT.TFile or pathlib.PurePath): Name of the root file
127  to open or an already open TFile instance
128  name (str): Name of the object in the file
129  title (str): Title to use for the object
130  contact (str): Contact person, usually in the form "Name <email>"
131  description (str): Text description what can be seen in the plot
132  check (str): Text description what to look for in the validation for
133  shifters to easily see if the distribution looks ok
134  xlabel (str): If given try to set this as the label for the x axis
135  ylabel (str): If given try to set this as the label for the y axis
136  metaoptions (str): Metaoptions (additional options to influence the
137  comparison between revisions, styling of the plot, etc.)
138 
139  .. warning::
140 
141  Different ways to specify LaTeX for different arguments:
142  see `create_validation_histograms`
143  """
144 
145  opened = False
146  if not isinstance(rootfile, ROOT.TFile):
147  if isinstance(rootfile, pathlib.PurePath):
148  rootfile = str(rootfile)
149  rootfile = ROOT.TFile(rootfile, "UPDATE")
150  opened = True
151  if not rootfile.IsOpen() or not rootfile.IsWritable():
152  raise RuntimeError(
153  f"ROOT file {rootfile.GetName()} is not open for writing"
154  )
155  obj = rootfile.Get(name)
156  if not obj:
157  raise RuntimeError(
158  f"Cannot find object named {name} in {rootfile.GetName()}"
159  )
160  validation_metadata_set(obj, *args, **argk)
161  # scope guard to avoid side effects by changing the global gDirectory
162  # in modules ...
163  # noinspection PyUnusedLocal
164  directory_guard = ROOT.TDirectory.TContext(rootfile)
165  obj.Write("", ROOT.TObject.kOverwrite)
166  if opened:
167  rootfile.Close()
168 
169 
170 class ValidationMetadataSetter(basf2.Module):
171  """
172  Simple module to set the valdiation metadata for a given list of objects
173  automatically at the end of event processing
174 
175  Just add this module **before** any
176  VariablesToNtuple/VariablesToHistogram modules and it will set the
177  correct validation metadata at the end of processing
178 
179  Warning:
180  The module needs to be before the modules creating the objects
181  as terminate() functions are executed in reverse order from last to
182  first module. If this module is after the creation modules the metadata
183  might not be set correctly
184  """
185 
186  def __init__(self, variables: List[Tuple[str]],
187  rootfile: Union[str, pathlib.PurePath], description=""):
188  """
189  Initialize ValidationMetadataSetter
190 
191  Arguments:
192  variables (list(tuple(str))): List of objects to set the metadata
193  for. Each entry should be the name of an object followed by the
194  metadata values which will be forwarded to
195  `validation_metadata_set`:
196  ``(name, title, contact, description, check, xlabel, ylabel,
197  metaoptions)``
198  where ``xlabel``, ``ylabel`` and ``metaoptions`` are optional
199  rootfile (str or pathlib.PurePath): The name of the ROOT file where
200  the objects can be found
201  description (str): Common description/information of/about all plots
202  in this ROOT file (will be displayed above the plots)
203  """
204  super().__init__()
205 
206  self._variables = variables
207  if isinstance(rootfile, pathlib.PurePath):
208  rootfile = str(rootfile)
209 
210  self._rootfile = rootfile # type: str
211 
213  self._description = description
214 
216  self._tfile = None # type: ROOT.TFile
217 
218  def initialize(self):
219  """Make sure we keep the file open"""
220  self._tfile = \
221  Belle2.RootFileCreationManager.getInstance().getFile(self._rootfile)
222 
223  def terminate(self):
224  """And update the metadata at the end"""
225  for name, *metadata in self._variables:
226  name = Belle2.makeROOTCompatible(name)
227  validation_metadata_update(self._tfile, name, *metadata)
228  if self._description:
229  file_description_set(self._tfile, self._description)
230  del self._tfile
231 
232 
233 def create_validation_histograms(
234  path: basf2.Path, rootfile: Union[str, pathlib.PurePath],
235  particlelist: str,
236  variables_1d: Optional[List[Tuple]] = None,
237  variables_2d: Optional[List[Tuple]] = None,
238  description=""
239 ) -> None:
240  """
241  Create histograms for all the variables and also label them to be useful
242  in validation plots in one go. This is similar to the
243  `modularAnalysis.variablesToHistogram` function but also sets the
244  metadata correctly to be used by the validation
245 
246  Arguments:
247  path (basf2.Path): Path where to put the modules
248  rootfile (str or pathlib.PurePath): Name of the output root file
249  particlelist (str): Name of the particle list, can be empty for event
250  dependent variables
251  variables_1d: List of 1D histogram definitions of the form
252  ``var, bins, min, max, title, contact, description, check_for
253  [, xlabel [, ylabel [, metaoptions]]]``
254  variables_2d: List of 2D histogram definitions of the form
255  ``var1, bins1, min1, max1, var2, bins2, min2, max2, title, contact,
256  description, check_for [, xlabel [, ylabel [, metaoptions]]]``
257  description: Common description/information of/about all plots in this
258  ROOT file (will be displayed above the plots)
259 
260  .. warning::
261 
262  Sadly, there are two different ways to specify latex formulas.
263 
264  1. The ROOT-style using ``#``:
265  ``"Einstein-Pythagoras a^{2} + b^{2} = #frac{E}{m}"``
266 
267  This style should be used for histogram title and labels, that is the
268  ``title``, ``xlabel`` and ``ylabel`` arguments
269 
270  2. The normal Latex style (with escaped backslashes or raw
271  string literals):
272  ``"Einstein-Pythagoras $a^2 + b^2 = \\\\frac{E}{m}$"``.
273 
274  This style should be used for all other fields like ``description``,
275  ``check_for``
276 
277  You can use the normal Latex style also for histogram title and descriptions
278  but the PDF version of the plots will still be buggy and not show all
279  formulas correctly.
280  """
281  if isinstance(rootfile, pathlib.PurePath):
282  rootfile = str(rootfile)
283  histograms_1d = []
284  histograms_2d = []
285  metadata = []
286  if variables_1d is not None:
287  for var, nbins, vmin, vmax, *data in variables_1d:
288  histograms_1d.append((var, nbins, vmin, vmax))
289  metadata.append([var] + data)
290 
291  if variables_2d is not None:
292  for row in variables_2d:
293  var1 = row[0]
294  var2 = row[4]
295  histograms_2d.append(row[:8])
296  metadata.append([var1 + var2] + list(row[8:]))
297 
298  path.add_module(
299  ValidationMetadataSetter(metadata, rootfile, description=description)
300  )
301  path.add_module(
302  "VariablesToHistogram",
303  particleList=particlelist,
304  variables=histograms_1d,
305  variables_2d=histograms_2d,
306  fileName=rootfile,
307  )
308 
309 # End suppression of doxygen checks
310 # @endcond
Belle2::makeROOTCompatible
std::string makeROOTCompatible(std::string str)
Remove special characters that ROOT dislikes in branch names, e.g.
Definition: MakeROOTCompatible.cc:74
Belle2::RootFileCreationManager::getInstance
static RootFileCreationManager & getInstance()
Interface for the FileManager.
Definition: RootFileCreationManager.cc:57