Logo ROOT   master
Reference Guide
InterfaceUtils.hxx
Go to the documentation of this file.
1 // Author: Enrico Guiraud, Danilo Piparo CERN 02/2018
2 
3 /*************************************************************************
4  * Copyright (C) 1995-2018, Rene Brun and Fons Rademakers. *
5  * All rights reserved. *
6  * *
7  * For the licensing terms see $ROOTSYS/LICENSE. *
8  * For the list of contributors see $ROOTSYS/README/CREDITS. *
9  *************************************************************************/
10 
11 #ifndef ROOT_RDF_TINTERFACE_UTILS
12 #define ROOT_RDF_TINTERFACE_UTILS
13 
14 #include <ROOT/RDF/RAction.hxx>
15 #include <ROOT/RDF/ActionHelpers.hxx> // for BuildAction
18 #include <ROOT/RDF/RFilter.hxx>
19 #include <ROOT/RDF/Utils.hxx>
25 #include <ROOT/RMakeUnique.hxx>
26 #include <ROOT/RStringView.hxx>
27 #include <ROOT/TypeTraits.hxx>
28 #include <TError.h> // gErrorIgnoreLevel
29 #include <TH1.h>
30 
31 #include <deque>
32 #include <functional>
33 #include <map>
34 #include <memory>
35 #include <string>
36 #include <type_traits>
37 #include <typeinfo>
38 #include <vector>
39 
40 class TObjArray;
41 class TTree;
42 namespace ROOT {
43 namespace Detail {
44 namespace RDF {
45 class RNodeBase;
46 }
47 }
48 namespace RDF {
49 template <typename T>
50 class RResultPtr;
51 template<typename T, typename V>
52 class RInterface;
54 class RDataSource;
55 } // namespace RDF
56 
57 } // namespace ROOT
58 
59 /// \cond HIDDEN_SYMBOLS
60 
61 namespace ROOT {
62 namespace Internal {
63 namespace RDF {
64 using namespace ROOT::Detail::RDF;
65 using namespace ROOT::RDF;
66 namespace TTraits = ROOT::TypeTraits;
68 
70 HeadNode_t CreateSnapshotRDF(const ColumnNames_t &validCols,
71  std::string_view treeName,
72  std::string_view fileName,
73  bool isLazy,
74  RLoopManager &loopManager,
75  std::unique_ptr<RDFInternal::RActionBase> actionPtr);
76 
77 std::string DemangleTypeIdName(const std::type_info &typeInfo);
78 
80  ROOT::RDF::RDataSource *dataSource, std::string_view columnNameRegexp,
81  std::string_view callerName);
82 
83 /// An helper object that sets and resets gErrorIgnoreLevel via RAII.
84 class RIgnoreErrorLevelRAII {
85 private:
86  int fCurIgnoreErrorLevel = gErrorIgnoreLevel;
87 
88 public:
89  RIgnoreErrorLevelRAII(int errorIgnoreLevel) { gErrorIgnoreLevel = errorIgnoreLevel; }
90  RIgnoreErrorLevelRAII() { gErrorIgnoreLevel = fCurIgnoreErrorLevel; }
91 };
92 
93 /****** BuildAction overloads *******/
94 
95 // clang-format off
96 /// This namespace defines types to be used for tag dispatching in RInterface.
97 namespace ActionTags {
98 struct Histo1D{};
99 struct Histo2D{};
100 struct Histo3D{};
101 struct Graph{};
102 struct Profile1D{};
103 struct Profile2D{};
104 struct Min{};
105 struct Max{};
106 struct Sum{};
107 struct Mean{};
108 struct Fill{};
109 struct StdDev{};
110 struct Display{};
111 }
112 // clang-format on
113 
114 template <typename T, bool ISV6HISTO = std::is_base_of<TH1, T>::value>
115 struct HistoUtils {
116  static void SetCanExtendAllAxes(T &h) { h.SetCanExtend(::TH1::kAllAxes); }
117  static bool HasAxisLimits(T &h)
118  {
119  auto xaxis = h.GetXaxis();
120  return !(xaxis->GetXmin() == 0. && xaxis->GetXmax() == 0.);
121  }
122 };
123 
124 template <typename T>
125 struct HistoUtils<T, false> {
126  static void SetCanExtendAllAxes(T &) {}
127  static bool HasAxisLimits(T &) { return true; }
128 };
129 
130 // Generic filling (covers Histo2D, Histo3D, Profile1D and Profile2D actions, with and without weights)
131 template <typename... BranchTypes, typename ActionTag, typename ActionResultType, typename PrevNodeType>
132 std::unique_ptr<RActionBase>
133 BuildAction(const ColumnNames_t &bl, const std::shared_ptr<ActionResultType> &h, const unsigned int nSlots,
134  std::shared_ptr<PrevNodeType> prevNode, ActionTag, RDFInternal::RBookedCustomColumns &&customColumns)
135 {
136  using Helper_t = FillParHelper<ActionResultType>;
137  using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<BranchTypes...>>;
138  return std::make_unique<Action_t>(Helper_t(h, nSlots), bl, std::move(prevNode), std::move(customColumns));
139 }
140 
141 // Histo1D filling (must handle the special case of distinguishing FillParHelper and FillHelper
142 template <typename... BranchTypes, typename PrevNodeType>
143 std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<::TH1D> &h,
144  const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
145  ActionTags::Histo1D, RDFInternal::RBookedCustomColumns &&customColumns)
146 {
147  auto hasAxisLimits = HistoUtils<::TH1D>::HasAxisLimits(*h);
148 
149  if (hasAxisLimits) {
150  using Helper_t = FillParHelper<::TH1D>;
151  using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<BranchTypes...>>;
152  return std::make_unique<Action_t>(Helper_t(h, nSlots), bl, std::move(prevNode), std::move(customColumns));
153  } else {
154  using Helper_t = FillHelper;
155  using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<BranchTypes...>>;
156  return std::make_unique<Action_t>(Helper_t(h, nSlots), bl, std::move(prevNode), std::move(customColumns));
157  }
158 }
159 
160 template <typename... BranchTypes, typename PrevNodeType>
161 std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<TGraph> &g,
162  const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
163  ActionTags::Graph, RDFInternal::RBookedCustomColumns &&customColumns)
164 {
165  using Helper_t = FillTGraphHelper;
166  using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<BranchTypes...>>;
167  return std::make_unique<Action_t>(Helper_t(g, nSlots), bl, std::move(prevNode), std::move(customColumns));
168 }
169 
170 // Min action
171 template <typename BranchType, typename PrevNodeType, typename ActionResultType>
172 std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<ActionResultType> &minV,
173  const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
175 {
176  using Helper_t = MinHelper<ActionResultType>;
177  using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<BranchType>>;
178  return std::make_unique<Action_t>(Helper_t(minV, nSlots), bl, std::move(prevNode), std::move(customColumns));
179 }
180 
181 // Max action
182 template <typename BranchType, typename PrevNodeType, typename ActionResultType>
183 std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<ActionResultType> &maxV,
184  const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
186 {
187  using Helper_t = MaxHelper<ActionResultType>;
188  using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<BranchType>>;
189  return std::make_unique<Action_t>(Helper_t(maxV, nSlots), bl, std::move(prevNode), std::move(customColumns));
190 }
191 
192 // Sum action
193 template <typename BranchType, typename PrevNodeType, typename ActionResultType>
194 std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<ActionResultType> &sumV,
195  const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
197 {
198  using Helper_t = SumHelper<ActionResultType>;
199  using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<BranchType>>;
200  return std::make_unique<Action_t>(Helper_t(sumV, nSlots), bl, std::move(prevNode), std::move(customColumns));
201 }
202 
203 // Mean action
204 template <typename BranchType, typename PrevNodeType>
205 std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<double> &meanV,
206  const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
208 {
209  using Helper_t = MeanHelper;
210  using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<BranchType>>;
211  return std::make_unique<Action_t>(Helper_t(meanV, nSlots), bl, std::move(prevNode), std::move(customColumns));
212 }
213 
214 // Standard Deviation action
215 template <typename BranchType, typename PrevNodeType>
216 std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<double> &stdDeviationV,
217  const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
219 {
220  using Helper_t = StdDevHelper;
221  using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<BranchType>>;
222  return std::make_unique<Action_t>(Helper_t(stdDeviationV, nSlots), bl, prevNode, std::move(customColumns));
223 }
224 
225 // Display action
226 template <typename... BranchTypes, typename PrevNodeType>
227 std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<RDisplay> &d,
228  const unsigned int, std::shared_ptr<PrevNodeType> prevNode,
229  ActionTags::Display, RDFInternal::RBookedCustomColumns &&customColumns)
230 {
231  using Helper_t = DisplayHelper<PrevNodeType>;
232  using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<BranchTypes...>>;
233  return std::make_unique<Action_t>(Helper_t(d, prevNode), bl, prevNode, std::move(customColumns));
234 }
235 
236 /****** end BuildAndBook ******/
237 
238 template <typename Filter>
239 void CheckFilter(Filter &)
240 {
241  using FilterRet_t = typename RDF::CallableTraits<Filter>::ret_type;
242  static_assert(std::is_convertible<FilterRet_t, bool>::value,
243  "filter expression returns a type that is not convertible to bool");
244 }
245 
246 void CheckCustomColumn(std::string_view definedCol, TTree *treePtr, const ColumnNames_t &customCols,
247  const std::map<std::string, std::string> &aliasMap, const ColumnNames_t &dataSourceColumns);
248 
249 std::string PrettyPrintAddr(const void *const addr);
250 
251 void BookFilterJit(RJittedFilter *jittedFilter, void *prevNodeOnHeap, std::string_view name,
252  std::string_view expression, const std::map<std::string, std::string> &aliasMap,
253  const ColumnNames_t &branches, const RDFInternal::RBookedCustomColumns &customCols, TTree *tree,
254  RDataSource *ds, unsigned int namespaceID);
255 
256 void BookDefineJit(std::string_view name, std::string_view expression, RLoopManager &lm, RDataSource *ds,
257  const std::shared_ptr<RJittedCustomColumn> &jittedCustomColumn,
258  const RDFInternal::RBookedCustomColumns &customCols, const ColumnNames_t &branches);
259 
260 std::string JitBuildAction(const ColumnNames_t &bl, void *prevNode, const std::type_info &art, const std::type_info &at,
261  void *r, TTree *tree, const unsigned int nSlots,
262  const RDFInternal::RBookedCustomColumns &customColumns, RDataSource *ds,
263  std::shared_ptr<RJittedAction> *jittedActionOnHeap, unsigned int namespaceID);
264 
265 // allocate a shared_ptr on the heap, return a reference to it. the user is responsible of deleting the shared_ptr*.
266 // this function is meant to only be used by RInterface's action methods, and should be deprecated as soon as we find
267 // a better way to make jitting work: the problem it solves is that we need to pass the same shared_ptr to the Helper
268 // object of each action and to the RResultPtr returned by the action. While the former is only instantiated when
269 // the event loop is about to start, the latter has to be returned to the user as soon as the action is booked.
270 // a heap allocated shared_ptr will stay alive long enough that at jitting time its address is still valid.
271 template <typename T>
272 std::shared_ptr<T> *MakeSharedOnHeap(const std::shared_ptr<T> &shPtr)
273 {
274  return new std::shared_ptr<T>(shPtr);
275 }
276 
277 bool AtLeastOneEmptyString(const std::vector<std::string_view> strings);
278 
279 /// Take a shared_ptr<AnyNodeType> and return a shared_ptr<RNodeBase>.
280 /// This works for RLoopManager nodes as well as filters and ranges.
281 std::shared_ptr<RNodeBase> UpcastNode(std::shared_ptr<RNodeBase> ptr);
282 
283 ColumnNames_t GetValidatedColumnNames(RLoopManager &lm, const unsigned int nColumns, const ColumnNames_t &columns,
284  const ColumnNames_t &validCustomColumns, RDataSource *ds);
285 
286 std::vector<bool> FindUndefinedDSColumns(const ColumnNames_t &requestedCols, const ColumnNames_t &definedDSCols);
287 
289 
290 template <typename T>
291 void AddDSColumnsHelper(RLoopManager &lm, std::string_view name, RDFInternal::RBookedCustomColumns &currentCols,
292  RDataSource &ds, unsigned int nSlots)
293 {
294  auto readers = ds.GetColumnReaders<T>(name);
295  auto getValue = [readers](unsigned int slot) { return *readers[slot]; };
297 
298  auto newCol = std::make_shared<NewCol_t>(&lm, name, std::move(getValue), ColumnNames_t{}, nSlots, currentCols,
299  /*isDSColumn=*/true);
300 
301  lm.RegisterCustomColumn(newCol.get());
302  currentCols.AddName(name);
303  currentCols.AddColumn(newCol, name);
304 }
305 
306 /// Take list of column names that must be defined, current map of custom columns, current list of defined column names,
307 /// and return a new map of custom columns (with the new datasource columns added to it)
308 template <typename... ColumnTypes, std::size_t... S>
310 AddDSColumns(RLoopManager &lm, const std::vector<std::string> &requiredCols,
311  const RDFInternal::RBookedCustomColumns &currentCols, RDataSource &ds, unsigned int nSlots,
312  std::index_sequence<S...>, TTraits::TypeList<ColumnTypes...>)
313 {
314 
315  const auto mustBeDefined = FindUndefinedDSColumns(requiredCols, currentCols.GetNames());
316  if (std::none_of(mustBeDefined.begin(), mustBeDefined.end(), [](bool b) { return b; })) {
317  // no need to define any column
318  return currentCols;
319  } else {
320  auto newColumns(currentCols);
321 
322  // hack to expand a template parameter pack without c++17 fold expressions.
323  int expander[] = {(mustBeDefined[S] ? AddDSColumnsHelper<ColumnTypes>(lm, requiredCols[S], newColumns, ds, nSlots)
324  : /*no-op*/ ((void)0),
325  0)...,
326  0};
327  (void)expander; // avoid unused variable warnings
328  (void)nSlots; // avoid unused variable warnings
329  return newColumns;
330  }
331 }
332 
333 // this function is meant to be called by the jitted code generated by BookFilterJit
334 template <typename F, typename PrevNode>
335 void JitFilterHelper(F &&f, const ColumnNames_t &cols, std::string_view name, RJittedFilter *jittedFilter,
336  std::shared_ptr<PrevNode> *prevNodeOnHeap, RDFInternal::RBookedCustomColumns *customColumns)
337 {
338  // mock Filter logic -- validity checks and Define-ition of RDataSource columns
339  using F_t = RFilter<F, PrevNode>;
340  using ColTypes_t = typename TTraits::CallableTraits<F>::arg_types;
341  constexpr auto nColumns = ColTypes_t::list_size;
342  RDFInternal::CheckFilter(f);
343 
344  auto &lm = *jittedFilter->GetLoopManagerUnchecked(); // RLoopManager must exist at this time
345  auto ds = lm.GetDataSource();
346 
347  auto newColumns = ds ? RDFInternal::AddDSColumns(lm, cols, *customColumns, *ds, lm.GetNSlots(),
348  std::make_index_sequence<nColumns>(), ColTypes_t())
349  : *customColumns;
350 
351  // customColumns points to the columns structure in the heap, created before the jitted call so that the jitter can
352  // share data after it has lazily compiled the code. Here the data has been used and the memory can be freed.
353  delete customColumns;
354 
355  jittedFilter->SetFilter(std::make_unique<F_t>(std::move(f), cols, *prevNodeOnHeap, newColumns, name));
356  delete prevNodeOnHeap;
357 }
358 
359 template <typename F>
360 void JitDefineHelper(F &&f, const ColumnNames_t &cols, std::string_view name, RLoopManager *lm,
361  RJittedCustomColumn &jittedCustomCol, RDFInternal::RBookedCustomColumns *customColumns)
362 {
364  using ColTypes_t = typename TTraits::CallableTraits<F>::arg_types;
365  constexpr auto nColumns = ColTypes_t::list_size;
366 
367  auto ds = lm->GetDataSource();
368  auto newColumns = ds ? RDFInternal::AddDSColumns(*lm, cols, *customColumns, *ds, lm->GetNSlots(),
369  std::make_index_sequence<nColumns>(), ColTypes_t())
370  : *customColumns;
371 
372  // customColumns points to the columns structure in the heap, created before the jitted call so that the jitter can
373  // share data after it has lazily compiled the code. Here the data has been used and the memory can be freed.
374  delete customColumns;
375 
376  // use unique_ptr<RCustomColumnBase> instead of make_unique<NewCol_t> to reduce jit/compile-times
377  jittedCustomCol.SetCustomColumn(
378  std::unique_ptr<RCustomColumnBase>(new NewCol_t(lm, name, std::move(f), cols, lm->GetNSlots(), newColumns)));
379 }
380 
381 /// Convenience function invoked by jitted code to build action nodes at runtime
382 template <typename ActionTag, typename... BranchTypes, typename PrevNodeType, typename ActionResultType>
383 void CallBuildAction(std::shared_ptr<PrevNodeType> *prevNodeOnHeap, const ColumnNames_t &bl, const unsigned int nSlots,
384  const std::shared_ptr<ActionResultType> *rOnHeap,
385  std::shared_ptr<RJittedAction> *jittedActionOnHeap,
386  RDFInternal::RBookedCustomColumns *customColumns)
387 {
388  // if we are here it means we are jitting, if we are jitting the loop manager must be alive
389  auto &prevNodePtr = *prevNodeOnHeap;
390  auto &loopManager = *prevNodePtr->GetLoopManagerUnchecked();
391  using ColTypes_t = TypeList<BranchTypes...>;
392  constexpr auto nColumns = ColTypes_t::list_size;
393  auto ds = loopManager.GetDataSource();
394  auto newColumns = ds ? RDFInternal::AddDSColumns(loopManager, bl, *customColumns, *ds, loopManager.GetNSlots(),
395  std::make_index_sequence<nColumns>(), ColTypes_t())
396  : *customColumns;
397 
398  auto actionPtr =
399  BuildAction<BranchTypes...>(bl, *rOnHeap, nSlots, std::move(prevNodePtr), ActionTag{}, std::move(newColumns));
400  (*jittedActionOnHeap)->SetAction(std::move(actionPtr));
401 
402  // customColumns points to the columns structure in the heap, created before the jitted call so that the jitter can
403  // share data after it has lazily compiled the code. Here the data has been used and the memory can be freed.
404  delete customColumns;
405 
406  delete rOnHeap;
407  delete prevNodeOnHeap;
408  delete jittedActionOnHeap;
409 }
410 
411 /// The contained `type` alias is `double` if `T == RInferredType`, `U` if `T == std::container<U>`, `T` otherwise.
412 template <typename T, bool Container = RDFInternal::IsDataContainer<T>::value && !std::is_same<T, std::string>::value>
413 struct TMinReturnType {
414  using type = T;
415 };
416 
417 template <>
418 struct TMinReturnType<RInferredType, false> {
419  using type = double;
420 };
421 
422 template <typename T>
423 struct TMinReturnType<T, true> {
425 };
426 
427 // return wrapper around f that prepends an `unsigned int slot` parameter
428 template <typename R, typename F, typename... Args>
429 std::function<R(unsigned int, Args...)> AddSlotParameter(F &f, TypeList<Args...>)
430 {
431  return [f](unsigned int, Args... a) -> R { return f(a...); };
432 }
433 
434 template <typename BranchType, typename... Rest>
435 struct TNeedJitting {
436  static constexpr bool value = TNeedJitting<Rest...>::value;
437 };
438 
439 template <typename... Rest>
440 struct TNeedJitting<RInferredType, Rest...> {
441  static constexpr bool value = true;
442 };
443 
444 template <typename T>
445 struct TNeedJitting<T> {
446  static constexpr bool value = false;
447 };
448 
449 template <>
450 struct TNeedJitting<RInferredType> {
451  static constexpr bool value = true;
452 };
453 
455 
456 ///////////////////////////////////////////////////////////////////////////////
457 /// Check preconditions for RInterface::Aggregate:
458 /// - the aggregator callable must have signature `U(U,T)` or `void(U&,T)`.
459 /// - the merge callable must have signature `U(U,U)` or `void(std::vector<U>&)`
461  typename mergeArgsNoDecay_t = typename CallableTraits<Merge>::arg_types_nodecay,
462  typename mergeArgs_t = typename CallableTraits<Merge>::arg_types,
463  typename mergeRet_t = typename CallableTraits<Merge>::ret_type>
464 void CheckAggregate(TypeList<U, T>)
465 {
466  constexpr bool isAggregatorOk =
467  (std::is_same<R, decayedU>::value) || (std::is_same<R, void>::value && std::is_lvalue_reference<U>::value);
468  static_assert(isAggregatorOk, "aggregator function must have signature `U(U,T)` or `void(U&,T)`");
469  constexpr bool isMergeOk =
470  (std::is_same<TypeList<decayedU, decayedU>, mergeArgs_t>::value && std::is_same<decayedU, mergeRet_t>::value) ||
471  (std::is_same<TypeList<std::vector<decayedU> &>, mergeArgsNoDecay_t>::value &&
472  std::is_same<void, mergeRet_t>::value);
473  static_assert(isMergeOk, "merge function must have signature `U(U,U)` or `void(std::vector<U>&)`");
474 }
475 
476 ///////////////////////////////////////////////////////////////////////////////
477 /// This overload of CheckAggregate is called when the aggregator takes more than two arguments
478 template <typename R, typename T>
479 void CheckAggregate(T)
480 {
481  static_assert(sizeof(T) == 0, "aggregator function must take exactly two arguments");
482 }
483 
484 ///////////////////////////////////////////////////////////////////////////////
485 /// Check as many template parameters were passed as the number of column names, throw if this is not the case.
486 void CheckTypesAndPars(unsigned int nTemplateParams, unsigned int nColumnNames);
487 
488 /// Return local BranchNames or default BranchNames according to which one should be used
489 const ColumnNames_t SelectColumns(unsigned int nArgs, const ColumnNames_t &bl, const ColumnNames_t &defBl);
490 
491 /// Check whether column names refer to a valid branch of a TTree or have been `Define`d. Return invalid column names.
492 ColumnNames_t FindUnknownColumns(const ColumnNames_t &requiredCols, const ColumnNames_t &datasetColumns,
493  const ColumnNames_t &definedCols, const ColumnNames_t &dataSourceColumns);
494 
495 bool IsInternalColumn(std::string_view colName);
496 
497 /// Returns the list of Filters defined in the whole graph
498 std::vector<std::string> GetFilterNames(const std::shared_ptr<RLoopManager> &loopManager);
499 
500 /// Returns the list of Filters defined in the branch
501 template <typename NodeType>
502 std::vector<std::string> GetFilterNames(const std::shared_ptr<NodeType> &node)
503 {
504  std::vector<std::string> filterNames;
505  node->AddFilterName(filterNames);
506  return filterNames;
507 }
508 
509 // Check if a condition is true for all types
510 template <bool...>
511 struct TBoolPack;
512 
513 template <bool... bs>
514 using IsTrueForAllImpl_t = typename std::is_same<TBoolPack<bs..., true>, TBoolPack<true, bs...>>;
515 
516 template <bool... Conditions>
517 struct TEvalAnd {
518  static constexpr bool value = IsTrueForAllImpl_t<Conditions...>::value;
519 };
520 
521 // Check if a class is a specialisation of stl containers templates
522 // clang-format off
523 
524 template <typename>
525 struct IsList_t : std::false_type {};
526 
527 template <typename T>
528 struct IsList_t<std::list<T>> : std::true_type {};
529 
530 template <typename>
531 struct IsDeque_t : std::false_type {};
532 
533 template <typename T>
534 struct IsDeque_t<std::deque<T>> : std::true_type {};
535 // clang-format on
536 
537 } // namespace RDF
538 } // namespace Internal
539 
540 namespace Detail {
541 namespace RDF {
542 
543 /// The aliased type is `double` if `T == RInferredType`, `U` if `T == container<U>`, `T` otherwise.
544 template <typename T>
545 using MinReturnType_t = typename RDFInternal::TMinReturnType<T>::type;
546 
547 template <typename T>
548 using MaxReturnType_t = MinReturnType_t<T>;
549 
550 template <typename T>
551 using SumReturnType_t = MinReturnType_t<T>;
552 
553 } // namespace RDF
554 } // namespace Detail
555 } // namespace ROOT
556 
557 /// \endcond
558 
559 #endif
Smart pointer for the return type of actions.
An array of TObjects.
Definition: TObjArray.h:37
The head node of a RDF computation graph.
ColumnNames_t GetTopLevelBranchNames(TTree &t)
Get all the top-level branches names, including the ones of the friend trees.
bool AtLeastOneEmptyString(const std::vector< std::string_view > strings)
R__EXTERN Int_t gErrorIgnoreLevel
Definition: TError.h:105
Returns the available number of logical cores.
Definition: StringConv.hxx:21
bool IsInternalColumn(std::string_view colName)
typename TakeFirstParameter< T >::type TakeFirstParameter_t
Definition: TypeTraits.hxx:155
double T(double x)
Definition: ChebyshevPol.h:34
#define g(i)
Definition: RSha256.hxx:105
void BookDefineJit(std::string_view name, std::string_view expression, RLoopManager &lm, RDataSource *ds, const std::shared_ptr< RJittedCustomColumn > &jittedCustomColumn, const RDFInternal::RBookedCustomColumns &customCols, const ColumnNames_t &branches)
#define f(i)
Definition: RSha256.hxx:104
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
STL namespace.
#define R(a, b, c, d, e, f, g, h, i)
Definition: RSha256.hxx:110
HeadNode_t CreateSnapshotRDF(const ColumnNames_t &validCols, std::string_view treeName, std::string_view fileName, bool isLazy, RLoopManager &loopManager, std::unique_ptr< RDFInternal::RActionBase > actionPtr)
std::shared_ptr< RNodeBase > UpcastNode(std::shared_ptr< RNodeBase > ptr)
const ColumnNames_t SelectColumns(unsigned int nRequiredNames, const ColumnNames_t &names, const ColumnNames_t &defaultNames)
Choose between local column names or default column names, throw in case of errors.
void SetCustomColumn(std::unique_ptr< RCustomColumnBase > c)
unsigned int GetNSlots() const
void function(const Char_t *name_, T fun, const Char_t *docstring=0)
Definition: RExports.h:151
A wrapper around a concrete RFilter, which forwards all calls to it RJittedFilter is the type of the ...
std::vector< std::string > GetFilterNames(const std::shared_ptr< RLoopManager > &loopManager)
RVec< T > Filter(const RVec< T > &v, F &&f)
Create a new collection with the elements passing the filter expressed by the predicate.
Definition: RVec.hxx:938
ColumnNames_t FindUnknownColumns(const ColumnNames_t &requiredCols, const ColumnNames_t &datasetColumns, const ColumnNames_t &definedCols, const ColumnNames_t &dataSourceColumns)
std::vector< bool > FindUndefinedDSColumns(const ColumnNames_t &requestedCols, const ColumnNames_t &definedCols)
Return a bitset each element of which indicates whether the corresponding element in selectedColumns ...
std::string DemangleTypeIdName(const std::type_info &typeInfo)
RooArgSet S(const RooAbsArg &v1)
#define F(x, y, z)
ROOT::R::TRInterface & r
Definition: Object.C:4
auto * a
Definition: textangle.C:12
void CheckCustomColumn(std::string_view definedCol, TTree *treePtr, const ColumnNames_t &customCols, const std::map< std::string, std::string > &aliasMap, const ColumnNames_t &dataSourceColumns)
void BookFilterJit(RJittedFilter *jittedFilter, void *prevNodeOnHeap, std::string_view name, std::string_view expression, const std::map< std::string, std::string > &aliasMap, const ColumnNames_t &branches, const RDFInternal::RBookedCustomColumns &customCols, TTree *tree, RDataSource *ds, unsigned int namespaceID)
Double_t Mean(Long64_t n, const T *a, const Double_t *w=0)
Return the weighted mean of an array a with length n.
Definition: TMath.h:1063
The public interface to the RDataFrame federation of classes.
#define h(i)
Definition: RSha256.hxx:106
Lightweight storage for a collection of types.
Definition: TypeTraits.hxx:25
Double_t StdDev(Long64_t n, const T *a, const Double_t *w=0)
Definition: TMath.h:518
virtual RLoopManager * GetLoopManagerUnchecked()
Definition: RNodeBase.hxx:64
#define d(i)
Definition: RSha256.hxx:102
ColumnNames_t GetValidatedColumnNames(RLoopManager &lm, const unsigned int nColumns, const ColumnNames_t &columns, const ColumnNames_t &validCustomColumns, RDataSource *ds)
Given the desired number of columns and the user-provided list of columns:
int type
Definition: TGX11.cxx:120
ROOT type_traits extensions.
Definition: TypeTraits.hxx:21
A wrapper around a concrete RCustomColumn, which forwards all calls to it RJittedCustomColumn is a pl...
Extract types from the signature of a callable object. See CallableTraits.
Definition: TypeTraits.hxx:36
T Sum(const RVec< T > &v)
Sum elements of an RVec.
Definition: RVec.hxx:758
typedef void((*Func_t)())
RDataSource * GetDataSource() const
std::string PrettyPrintAddr(const void *const addr)
void CheckTypesAndPars(unsigned int nTemplateParams, unsigned int nColumnNames)
void SetFilter(std::unique_ptr< RFilterBase > f)
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
Definition: TRolke.cxx:630
std::vector< T ** > GetColumnReaders(std::string_view columnName)
Called at most once per column by RDF.
Definition: tree.py:1
Encapsulates the columns defined by the user.
A TTree represents a columnar dataset.
Definition: TTree.h:72
RDataSource defines an API that RDataFrame can use to read arbitrary data formats.
std::string JitBuildAction(const ColumnNames_t &bl, void *prevNode, const std::type_info &art, const std::type_info &at, void *rOnHeap, TTree *tree, const unsigned int nSlots, const RDFInternal::RBookedCustomColumns &customCols, RDataSource *ds, std::shared_ptr< RJittedAction > *jittedActionOnHeap, unsigned int namespaceID)
ColumnNames_t GetNames() const
Returns the list of the names of the defined columns.
char name[80]
Definition: TGX11.cxx:109
ROOT::Detail::RDF::ColumnNames_t ColumnNames_t
Definition: RDataFrame.cxx:797
ColumnNames_t ConvertRegexToColumns(const RDFInternal::RBookedCustomColumns &customColumns, TTree *tree, ROOT::RDF::RDataSource *dataSource, std::string_view columnNameRegexp, std::string_view callerName)