Logo ROOT   6.14/05
Reference Guide
RDFNodes.hxx
Go to the documentation of this file.
1 // Author: Enrico Guiraud, Danilo Piparo CERN 03/2017
2 
3 /*************************************************************************
4  * Copyright (C) 1995-2016, 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_RDFNODES
12 #define ROOT_RDFNODES
13 
14 #ifndef NDEBUG
15 #include "TError.h"
16 #endif
18 #include "ROOT/TypeTraits.hxx"
19 #include "ROOT/RCutFlowReport.hxx"
20 #include "ROOT/RDataSource.hxx"
21 #include "ROOT/RDFNodesUtils.hxx"
22 #include "ROOT/RDFUtils.hxx"
23 #include "ROOT/RVec.hxx"
24 #include "ROOT/TSpinMutex.hxx"
25 #include "TTreeReaderArray.h"
26 #include "TTreeReaderValue.h"
27 #include "TError.h"
28 
29 #include <map>
30 #include <numeric> // std::accumulate (FillReport), std::iota (TSlotStack)
31 #include <string>
32 #include <tuple>
33 #include <cassert>
34 #include <climits>
35 #include <deque> // std::vector substitute in case of vector<bool>
36 #include <functional>
37 
38 namespace ROOT {
39 namespace Internal {
40 namespace RDF {
41 class RActionBase;
42 
43 // This is an helper class to allow to pick a slot without resorting to a map
44 // indexed by thread ids.
45 // WARNING: this class does not work as a regular stack. The size is
46 // fixed at construction time and no blocking is foreseen.
47 class TSlotStack {
48 private:
49  unsigned int &GetCount()
50  {
51  TTHREAD_TLS(unsigned int) count = 0U;
52  return count;
53  }
54  unsigned int &GetIndex()
55  {
56  TTHREAD_TLS(unsigned int) index = UINT_MAX;
57  return index;
58  }
59  unsigned int fCursor;
60  std::vector<unsigned int> fBuf;
62 
63 public:
64  TSlotStack() = delete;
65  TSlotStack(unsigned int size) : fCursor(size), fBuf(size) { std::iota(fBuf.begin(), fBuf.end(), 0U); }
66  void ReturnSlot(unsigned int slotNumber);
67  unsigned int GetSlot();
68 };
69 } // ns RDF
70 } // ns Internal
71 
72 namespace Detail {
73 namespace RDF {
74 using namespace ROOT::TypeTraits;
76 
77 // forward declarations for RLoopManager
78 using ActionBasePtr_t = std::shared_ptr<RDFInternal::RActionBase>;
79 using ActionBaseVec_t = std::vector<ActionBasePtr_t>;
80 class RCustomColumnBase;
81 using RCustomColumnBasePtr_t = std::shared_ptr<RCustomColumnBase>;
82 class RFilterBase;
83 using FilterBasePtr_t = std::shared_ptr<RFilterBase>;
84 using FilterBaseVec_t = std::vector<FilterBasePtr_t>;
85 class RRangeBase;
86 using RangeBasePtr_t = std::shared_ptr<RRangeBase>;
87 using RangeBaseVec_t = std::vector<RangeBasePtr_t>;
88 
89 class RLoopManager {
91  enum class ELoopType { kROOTFiles, kROOTFilesMT, kNoFiles, kNoFilesMT, kDataSource, kDataSourceMT };
92 
93  using Callback_t = std::function<void(unsigned int)>;
94  class TCallback {
97  std::vector<ULong64_t> fCounters;
98 
99  public:
100  TCallback(ULong64_t everyN, Callback_t &&f, unsigned int nSlots)
101  : fFun(std::move(f)), fEveryN(everyN), fCounters(nSlots, 0ull)
102  {
103  }
104 
105  void operator()(unsigned int slot)
106  {
107  auto &c = fCounters[slot];
108  ++c;
109  if (c == fEveryN) {
110  c = 0ull;
111  fFun(slot);
112  }
113  }
114  };
115 
118  std::vector<int> fHasBeenCalled; // std::vector<bool> is thread-unsafe for our purposes (and generally evil)
119 
120  public:
121  TOneTimeCallback(Callback_t &&f, unsigned int nSlots) : fFun(std::move(f)), fHasBeenCalled(nSlots, 0) {}
122 
123  void operator()(unsigned int slot)
124  {
125  if (fHasBeenCalled[slot] == 1)
126  return;
127  fFun(slot);
128  fHasBeenCalled[slot] = 1;
129  }
130  };
131 
134  FilterBaseVec_t fBookedNamedFilters; ///< Contains a subset of fBookedFilters, i.e. only the named filters
135  std::map<std::string, RCustomColumnBasePtr_t> fBookedCustomColumns;
136  ColumnNames_t fCustomColumnNames; ///< Contains names of all custom columns defined in the functional graph.
138  std::vector<std::shared_ptr<bool>> fResProxyReadiness;
139  ::TDirectory *const fDirPtr{nullptr};
140  std::shared_ptr<TTree> fTree{nullptr}; //< Shared pointer to the input TTree/TChain. It does not own the pointee if
141  // the TTree/TChain was passed directly as an argument to RDataFrame's ctor (in
142  // which case we let users retain ownership).
143  const ColumnNames_t fDefaultColumns;
144  const ULong64_t fNEmptyEntries{0};
145  const unsigned int fNSlots{1};
146  bool fMustRunNamedFilters{true};
147  unsigned int fNChildren{0}; ///< Number of nodes of the functional graph hanging from this object
148  unsigned int fNStopsReceived{0}; ///< Number of times that a children node signaled to stop processing entries.
149  const ELoopType fLoopType; ///< The kind of event loop that is going to be run (e.g. on ROOT files, on no files)
150  std::string fToJit; ///< string containing all `BuildAndBook` actions that should be jitted before running
151  const std::unique_ptr<RDataSource> fDataSource; ///< Owning pointer to a data-source object. Null if no data-source
152  ColumnNames_t fDefinedDataSourceColumns; ///< List of data-source columns that have been `Define`d so far
153  std::map<std::string, std::string> fAliasColumnNameMap; ///< ColumnNameAlias-columnName pairs
154  std::vector<TCallback> fCallbacks; ///< Registered callbacks
155  std::vector<TOneTimeCallback> fCallbacksOnce; ///< Registered callbacks to invoke just once before running the loop
156  /// A unique ID that identifies the computation graph that starts with this RLoopManager.
157  /// Used, for example, to jit objects in a namespace reserved for this computation graph
158  const unsigned int fID = GetNextID();
159 
160  void RunEmptySourceMT();
161  void RunEmptySource();
162  void RunTreeProcessorMT();
163  void RunTreeReader();
164  void RunDataSourceMT();
165  void RunDataSource();
166  void RunAndCheckFilters(unsigned int slot, Long64_t entry);
167  void InitNodeSlots(TTreeReader *r, unsigned int slot);
168  void InitNodes();
169  void CleanUpNodes();
170  void CleanUpTask(unsigned int slot);
171  void JitActions();
172  void EvalChildrenCounts();
173  unsigned int GetNextID() const;
174 
175 public:
176  RLoopManager(TTree *tree, const ColumnNames_t &defaultBranches);
177  RLoopManager(ULong64_t nEmptyEntries);
178  RLoopManager(std::unique_ptr<RDataSource> ds, const ColumnNames_t &defaultBranches);
179  RLoopManager(const RLoopManager &) = delete;
180  RLoopManager &operator=(const RLoopManager &) = delete;
181 
182  void Run();
183  RLoopManager *GetLoopManagerUnchecked();
184  const ColumnNames_t &GetDefaultColumnNames() const;
185  const ColumnNames_t &GetCustomColumnNames() const { return fCustomColumnNames; };
186  TTree *GetTree() const;
187  const std::map<std::string, RCustomColumnBasePtr_t> &GetBookedColumns() const { return fBookedCustomColumns; }
188  ::TDirectory *GetDirectory() const;
189  ULong64_t GetNEmptyEntries() const { return fNEmptyEntries; }
190  RDataSource *GetDataSource() const { return fDataSource.get(); }
191  void Book(const ActionBasePtr_t &actionPtr);
192  void Book(const FilterBasePtr_t &filterPtr);
193  void Book(const RCustomColumnBasePtr_t &branchPtr);
194  void Book(const std::shared_ptr<bool> &branchPtr);
195  void Book(const RangeBasePtr_t &rangePtr);
196  bool CheckFilters(int, unsigned int);
197  unsigned int GetNSlots() const { return fNSlots; }
198  bool MustRunNamedFilters() const { return fMustRunNamedFilters; }
199  void Report(ROOT::RDF::RCutFlowReport &rep) const;
200  /// End of recursive chain of calls, does nothing
202  void SetTree(const std::shared_ptr<TTree> &tree) { fTree = tree; }
203  void IncrChildrenCount() { ++fNChildren; }
204  void StopProcessing() { ++fNStopsReceived; }
205  void ToJit(const std::string &s) { fToJit.append(s); }
206  const ColumnNames_t &GetDefinedDataSourceColumns() const { return fDefinedDataSourceColumns; }
207  void AddDataSourceColumn(std::string_view name) { fDefinedDataSourceColumns.emplace_back(name); }
208  void AddColumnAlias(const std::string &alias, const std::string &colName) { fAliasColumnNameMap[alias] = colName; }
209  void AddCustomColumnName(std::string_view name) { fCustomColumnNames.emplace_back(name); }
210  const std::map<std::string, std::string> &GetAliasMap() const { return fAliasColumnNameMap; }
211  void RegisterCallback(ULong64_t everyNEvents, std::function<void(unsigned int)> &&f);
212  unsigned int GetID() const { return fID; }
213 };
214 } // end ns RDF
215 } // end ns Detail
216 
217 namespace Internal {
218 namespace RDF {
219 using namespace ROOT::Detail::RDF;
220 
221 /**
222 \class ROOT::Internal::RDF::TColumnValue
223 \ingroup dataframe
224 \brief Helper class that updates and returns TTree branches as well as RDataFrame temporary columns
225 \tparam T The type of the column
226 
227 RDataFrame nodes must access two different types of values during the event loop:
228 values of real branches, for which TTreeReader{Values,Arrays} act as proxies, or
229 temporary columns whose values are generated on the fly. While the type of the
230 value is known at compile time (or just-in-time), it is only at runtime that nodes
231 can check whether a certain value is generated on the fly or not.
232 
233 TColumnValue abstracts this difference by providing the same interface for
234 both cases and handling the reading or generation of new values transparently.
235 Only one of the two data members fReaderProxy or fValuePtr will be non-null
236 for a given TColumnValue, depending on whether the value comes from a real
237 TTree branch or from a temporary column respectively.
238 
239 RDataFrame nodes can store tuples of TColumnValues and retrieve an updated
240 value for the column via the `Get` method.
241 **/
242 template <typename T, bool MustUseRVec = IsRVec_t<T>::value>
244  // ColumnValue_t is the type of the column or the type of the elements of an array column
245  using ColumnValue_t = typename std::conditional<MustUseRVec, TakeFirstParameter_t<T>, T>::type;
246  using TreeReader_t =
247  typename std::conditional<MustUseRVec, TTreeReaderArray<ColumnValue_t>, TTreeReaderValue<ColumnValue_t>>::type;
248 
249  /// TColumnValue has a slightly different behaviour whether the column comes from a TTreeReader, a RDataFrame Define
250  /// or a RDataSource. It stores which it is as an enum.
251  enum class EColumnKind { kTree, kCustomColumn, kDataSource, kInvalid };
252  // Set to the correct value by MakeProxy or SetTmpColumn
254  /// The slot this value belongs to. Only needed when querying custom column values, it is set in `SetTmpColumn`.
255  unsigned int fSlot = std::numeric_limits<unsigned int>::max();
256 
257  // Each element of the following data members will be in use by a _single task_.
258  // The vectors are used as very small stacks (1-2 elements typically) that fill in case of interleaved task execution
259  // i.e. when more than one task needs readers in this worker thread.
260 
261  /// Owning ptrs to a TTreeReaderValue or TTreeReaderArray. Only used for Tree columns.
262  std::vector<std::unique_ptr<TreeReader_t>> fTreeReaders;
263  /// Non-owning ptrs to the value of a custom column.
264  std::vector<T *> fCustomValuePtrs;
265  /// Non-owning ptrs to the value of a data-source column.
266  std::vector<T **> fDSValuePtrs;
267  /// Non-owning ptrs to the node responsible for the custom column. Needed when querying custom values.
268  std::vector<RCustomColumnBase *> fCustomColumns;
269  /// Enumerator for the different properties of the branch storage in memory
270  enum class EStorageType : char { kContiguous, kUnknown, kSparse};
271  /// Signal whether we ever checked that the branch we are reading with a TTreeReaderArray stores array elements
272  /// in contiguous memory. Only used when T == RVec<U>.
274  /// If MustUseRVec, i.e. we are reading an array, we return a reference to this RVec to clients
276  bool fCopyWarningPrinted = false;
277 
278 public:
279  static constexpr bool fgMustUseRVec = MustUseRVec;
280 
281  TColumnValue() = default;
282 
283  void SetTmpColumn(unsigned int slot, RCustomColumnBase *tmpColumn);
284 
285  void MakeProxy(TTreeReader *r, const std::string &bn)
286  {
287  fColumnKind = EColumnKind::kTree;
288  fTreeReaders.emplace_back(new TreeReader_t(*r, bn.c_str()));
289  }
290 
291  /// This overload is used to return scalar quantities (i.e. types that are not read into a RVec)
292  template <typename U = T, typename std::enable_if<!TColumnValue<U>::fgMustUseRVec, int>::type = 0>
293  T &Get(Long64_t entry);
294 
295  /// This overload is used to return arrays (i.e. types that are read into a RVec).
296  /// In this case the returned T is always a RVec<ColumnValue_t>.
297  template <typename U = T, typename std::enable_if<TColumnValue<U>::fgMustUseRVec, int>::type = 0>
298  T &Get(Long64_t entry);
299 
300  void Reset()
301  {
302  switch (fColumnKind) {
303  case EColumnKind::kTree: fTreeReaders.pop_back(); break;
304  case EColumnKind::kCustomColumn:
305  fCustomColumns.pop_back();
306  fCustomValuePtrs.pop_back();
307  break;
308  case EColumnKind::kDataSource:
309  fCustomColumns.pop_back();
310  fDSValuePtrs.pop_back();
311  break;
312  case EColumnKind::kInvalid: throw std::runtime_error("ColumnKind not set for this TColumnValue");
313  }
314  }
315 };
316 
317 template <typename T>
319 };
320 
321 template <typename... BranchTypes>
322 struct TRDFValueTuple<TypeList<BranchTypes...>> {
323  using type = std::tuple<TColumnValue<BranchTypes>...>;
324 };
325 
326 template <typename BranchType>
328 
329 /// Clear the proxies of a tuple of TColumnValues
330 template <typename ValueTuple, std::size_t... S>
332 {
333  // hack to expand a parameter pack without c++17 fold expressions.
334  std::initializer_list<int> expander{(std::get<S>(values).Reset(), 0)...};
335  (void)expander; // avoid "unused variable" warnings
336 }
337 
338 class RActionBase {
339 protected:
340  RLoopManager *fLoopManager; ///< A raw pointer to the RLoopManager at the root of this functional
341  /// graph. It is only guaranteed to contain a valid address during an
342  /// event loop.
343  const unsigned int fNSlots; ///< Number of thread slots used by this node.
344 
345 public:
346  RActionBase(RLoopManager *implPtr, const unsigned int nSlots);
347  RActionBase(const RActionBase &) = delete;
348  RActionBase &operator=(const RActionBase &) = delete;
349  virtual ~RActionBase() = default;
350 
351  virtual void Run(unsigned int slot, Long64_t entry) = 0;
352  virtual void Initialize() = 0;
353  virtual void InitSlot(TTreeReader *r, unsigned int slot) = 0;
354  virtual void TriggerChildrenCount() = 0;
355  virtual void FinalizeSlot(unsigned int) = 0;
356  /// This method is invoked to update a partial result during the event loop, right before passing the result to a
357  /// user-defined callback registered via RResultPtr::RegisterCallback
358  virtual void *PartialUpdate(unsigned int slot) = 0;
359 };
360 
361 template <typename Helper, typename PrevDataFrame, typename ColumnTypes_t = typename Helper::ColumnTypes_t>
362 class RAction final : public RActionBase {
364 
365  Helper fHelper;
366  const ColumnNames_t fBranches;
367  PrevDataFrame &fPrevData;
368  std::vector<RDFValueTuple_t<ColumnTypes_t>> fValues;
369 
370 public:
371  RAction(Helper &&h, const ColumnNames_t &bl, PrevDataFrame &pd)
372  : RActionBase(pd.GetLoopManagerUnchecked(), pd.GetLoopManagerUnchecked()->GetNSlots()), fHelper(std::move(h)),
373  fBranches(bl), fPrevData(pd), fValues(fNSlots)
374  {
375  }
376 
377  RAction(const RAction &) = delete;
378  RAction &operator=(const RAction &) = delete;
379  ~RAction() { fHelper.Finalize(); }
380 
381  void Initialize() final { fHelper.Initialize(); }
382 
383  void InitSlot(TTreeReader *r, unsigned int slot) final
384  {
385  InitRDFValues(slot, fValues[slot], r, fBranches, fLoopManager->GetCustomColumnNames(),
386  fLoopManager->GetBookedColumns(), TypeInd_t());
387  fHelper.InitTask(r, slot);
388  }
389 
390  void Run(unsigned int slot, Long64_t entry) final
391  {
392  // check if entry passes all filters
393  if (fPrevData.CheckFilters(slot, entry))
394  Exec(slot, entry, TypeInd_t());
395  }
396 
397  template <std::size_t... S>
398  void Exec(unsigned int slot, Long64_t entry, std::index_sequence<S...>)
399  {
400  (void)entry; // avoid bogus 'unused parameter' warning in gcc4.9
401  fHelper.Exec(slot, std::get<S>(fValues[slot]).Get(entry)...);
402  }
403 
404  void TriggerChildrenCount() final { fPrevData.IncrChildrenCount(); }
405 
406  void FinalizeSlot(unsigned int slot) final
407  {
408  ClearValueReaders(slot);
409  fHelper.CallFinalizeTask(slot);
410  }
411 
412  void ClearValueReaders(unsigned int slot) { ResetRDFValueTuple(fValues[slot], TypeInd_t()); }
413 
414  /// This method is invoked to update a partial result during the event loop, right before passing the result to a
415  /// user-defined callback registered via RResultPtr::RegisterCallback
416  /// TODO the PartialUpdateImpl trick can go away once all action helpers will implement PartialUpdate
417  void *PartialUpdate(unsigned int slot) final { return PartialUpdateImpl(slot); }
418 
419 private:
420  // this overload is SFINAE'd out if Helper does not implement `PartialUpdate`
421  // the template parameter is required to defer instantiation of the method to SFINAE time
422  template <typename H = Helper>
423  auto PartialUpdateImpl(unsigned int slot) -> decltype(std::declval<H>().PartialUpdate(slot), (void *)(nullptr))
424  {
425  return &fHelper.PartialUpdate(slot);
426  }
427  // this one is always available but has lower precedence thanks to `...`
428  void *PartialUpdateImpl(...) { throw std::runtime_error("This action does not support callbacks yet!"); }
429 };
430 
431 } // end NS RDF
432 } // end NS Internal
433 
434 namespace Detail {
435 namespace RDF {
436 
438 protected:
439  RLoopManager *fLoopManager; ///< A raw pointer to the RLoopManager at the root of this functional graph. It is only
440  /// guaranteed to contain a valid address during an event loop.
441  const std::string fName;
442  unsigned int fNChildren{0}; ///< number of nodes of the functional graph hanging from this object
443  unsigned int fNStopsReceived{0}; ///< number of times that a children node signaled to stop processing entries.
444  const unsigned int fNSlots; ///< number of thread slots used by this node, inherited from parent node.
445  const bool fIsDataSourceColumn; ///< does the custom column refer to a data-source column? (or a user-define column?)
446  std::vector<Long64_t> fLastCheckedEntry;
447 
448 public:
449  RCustomColumnBase(RLoopManager *df, std::string_view name, const unsigned int nSlots, const bool isDSColumn);
450  RCustomColumnBase &operator=(const RCustomColumnBase &) = delete;
451  virtual ~RCustomColumnBase(); // outlined defaulted.
452  virtual void InitSlot(TTreeReader *r, unsigned int slot) = 0;
453  virtual void *GetValuePtr(unsigned int slot) = 0;
454  virtual const std::type_info &GetTypeId() const = 0;
455  RLoopManager *GetLoopManagerUnchecked() const;
456  std::string GetName() const;
457  virtual void Update(unsigned int slot, Long64_t entry) = 0;
458  virtual void ClearValueReaders(unsigned int slot) = 0;
459  bool IsDataSourceColumn() const { return fIsDataSourceColumn; }
460  void InitNode();
461 };
462 
463 // clang-format off
464 namespace TCCHelperTypes {
465 struct TNothing;
466 struct TSlot;
467 struct TSlotAndEntry;
468 }
469 // clang-format on
470 
471 template <typename F, typename UPDATE_HELPER_TYPE = TCCHelperTypes::TNothing>
472 class RCustomColumn final : public RCustomColumnBase {
473  // shortcuts
474  using TNothing = TCCHelperTypes::TNothing;
475  using TSlot = TCCHelperTypes::TSlot;
476  using TSlotAndEntry = TCCHelperTypes::TSlotAndEntry;
477  using UHT_t = UPDATE_HELPER_TYPE;
478  // other types
479  using FunParamTypes_t = typename CallableTraits<F>::arg_types;
480  using BranchTypesTmp_t =
481  typename RDFInternal::RemoveFirstParameterIf<std::is_same<TSlot, UHT_t>::value, FunParamTypes_t>::type;
482  using ColumnTypes_t = typename RDFInternal::RemoveFirstTwoParametersIf<std::is_same<TSlotAndEntry, UHT_t>::value,
485  using ret_type = typename CallableTraits<F>::ret_type;
486  // Avoid instantiating vector<bool> as `operator[]` returns temporaries in that case. Use std::deque instead.
487  using ValuesPerSlot_t =
488  typename std::conditional<std::is_same<ret_type, bool>::value, std::deque<ret_type>, std::vector<ret_type>>::type;
489 
491  const ColumnNames_t fBranches;
493 
494  std::vector<RDFInternal::RDFValueTuple_t<ColumnTypes_t>> fValues;
495 
496 public:
497  RCustomColumn(std::string_view name, F &&expression, const ColumnNames_t &bl, RLoopManager *lm,
498  bool isDSColumn = false)
499  : RCustomColumnBase(lm, name, lm->GetNSlots(), isDSColumn), fExpression(std::move(expression)), fBranches(bl),
500  fLastResults(fNSlots), fValues(fNSlots)
501  {
502  }
503 
504  RCustomColumn(const RCustomColumn &) = delete;
505  RCustomColumn &operator=(const RCustomColumn &) = delete;
506 
507  void InitSlot(TTreeReader *r, unsigned int slot) final
508  {
509  RDFInternal::InitRDFValues(slot, fValues[slot], r, fBranches, fLoopManager->GetCustomColumnNames(),
510  fLoopManager->GetBookedColumns(), TypeInd_t());
511  }
512 
513  void *GetValuePtr(unsigned int slot) final { return static_cast<void *>(&fLastResults[slot]); }
514 
515  void Update(unsigned int slot, Long64_t entry) final
516  {
517  if (entry != fLastCheckedEntry[slot]) {
518  // evaluate this filter, cache the result
519  UpdateHelper(slot, entry, TypeInd_t(), ColumnTypes_t(), (UPDATE_HELPER_TYPE *)nullptr);
520  fLastCheckedEntry[slot] = entry;
521  }
522  }
523 
524  const std::type_info &GetTypeId() const
525  {
526  return fIsDataSourceColumn ? typeid(typename std::remove_pointer<ret_type>::type) : typeid(ret_type);
527  }
528 
529  template <std::size_t... S, typename... BranchTypes>
531  TCCHelperTypes::TNothing *)
532  {
533  fLastResults[slot] = fExpression(std::get<S>(fValues[slot]).Get(entry)...);
534  // silence "unused parameter" warnings in gcc
535  (void)slot;
536  (void)entry;
537  }
538 
539  template <std::size_t... S, typename... BranchTypes>
541  TCCHelperTypes::TSlot *)
542  {
543  fLastResults[slot] = fExpression(slot, std::get<S>(fValues[slot]).Get(entry)...);
544  // silence "unused parameter" warnings in gcc
545  (void)slot;
546  (void)entry;
547  }
548 
549  template <std::size_t... S, typename... BranchTypes>
551  TCCHelperTypes::TSlotAndEntry *)
552  {
553  fLastResults[slot] = fExpression(slot, entry, std::get<S>(fValues[slot]).Get(entry)...);
554  // silence "unused parameter" warnings in gcc
555  (void)slot;
556  (void)entry;
557  }
558 
559  void ClearValueReaders(unsigned int slot) final { RDFInternal::ResetRDFValueTuple(fValues[slot], TypeInd_t()); }
560 };
561 
562 class RFilterBase {
563 protected:
564  RLoopManager *fLoopManager; ///< A raw pointer to the RLoopManager at the root of this functional graph. It is only
565  /// guaranteed to contain a valid address during an event loop.
566  std::vector<Long64_t> fLastCheckedEntry;
567  std::vector<int> fLastResult = {true}; // std::vector<bool> cannot be used in a MT context safely
568  std::vector<ULong64_t> fAccepted = {0};
569  std::vector<ULong64_t> fRejected = {0};
570  const std::string fName;
571  unsigned int fNChildren{0}; ///< Number of nodes of the functional graph hanging from this object
572  unsigned int fNStopsReceived{0}; ///< Number of times that a children node signaled to stop processing entries.
573  const unsigned int fNSlots; ///< Number of thread slots used by this node, inherited from parent node.
574 
575 public:
576  RFilterBase(RLoopManager *df, std::string_view name, const unsigned int nSlots);
577  RFilterBase &operator=(const RFilterBase &) = delete;
578  virtual ~RFilterBase() = default;
579 
580  virtual void InitSlot(TTreeReader *r, unsigned int slot) = 0;
581  virtual bool CheckFilters(unsigned int slot, Long64_t entry) = 0;
582  virtual void Report(ROOT::RDF::RCutFlowReport &) const = 0;
583  virtual void PartialReport(ROOT::RDF::RCutFlowReport &) const = 0;
584  RLoopManager *GetLoopManagerUnchecked() const;
585  bool HasName() const;
586  virtual void FillReport(ROOT::RDF::RCutFlowReport &) const;
587  virtual void IncrChildrenCount() = 0;
588  virtual void StopProcessing() = 0;
589  virtual void ResetChildrenCount()
590  {
591  fNChildren = 0;
592  fNStopsReceived = 0;
593  }
594  virtual void TriggerChildrenCount() = 0;
595  virtual void ResetReportCount()
596  {
597  assert(!fName.empty()); // this method is to only be called on named filters
598  // fAccepted and fRejected could be different than 0 if this is not the first event-loop run using this filter
599  std::fill(fAccepted.begin(), fAccepted.end(), 0);
600  std::fill(fRejected.begin(), fRejected.end(), 0);
601  }
602  virtual void ClearValueReaders(unsigned int slot) = 0;
603  virtual void InitNode();
604 };
605 
606 /// A wrapper around a concrete RFilter, which forwards all calls to it
607 /// RJittedFilter is the type of the node returned by jitted Filter calls: the concrete filter can be created and set
608 /// at a later time, from jitted code.
609 // FIXME after switching to the new ownership model, RJittedFilter should avoid inheriting from RFilterBase and
610 // overriding all of its methods: it can just implement them, and RFilterBase's can go back to have non-virtual methods
611 class RJittedFilter final : public RFilterBase {
612  std::unique_ptr<RFilterBase> fConcreteFilter = nullptr;
613 
614 public:
616 
617  void SetFilter(std::unique_ptr<RFilterBase> f);
618 
619  void InitSlot(TTreeReader *r, unsigned int slot) override final;
620  bool CheckFilters(unsigned int slot, Long64_t entry) override final;
621  void Report(ROOT::RDF::RCutFlowReport &) const override final;
622  void PartialReport(ROOT::RDF::RCutFlowReport &) const override final;
623  void FillReport(ROOT::RDF::RCutFlowReport &) const override final;
624  void IncrChildrenCount() override final;
625  void StopProcessing() override final;
626  void ResetChildrenCount() override final;
627  void TriggerChildrenCount() override final;
628  void ResetReportCount() override final;
629  void ClearValueReaders(unsigned int slot) override final;
630  void InitNode() override final;
631 };
632 
633 template <typename FilterF, typename PrevDataFrame>
634 class RFilter final : public RFilterBase {
635  using ColumnTypes_t = typename CallableTraits<FilterF>::arg_types;
637 
638  FilterF fFilter;
639  const ColumnNames_t fBranches;
640  PrevDataFrame &fPrevData;
641  std::vector<RDFInternal::RDFValueTuple_t<ColumnTypes_t>> fValues;
642 
643 public:
644  RFilter(FilterF &&f, const ColumnNames_t &bl, PrevDataFrame &pd, std::string_view name = "")
645  : RFilterBase(pd.GetLoopManagerUnchecked(), name, pd.GetLoopManagerUnchecked()->GetNSlots()),
646  fFilter(std::move(f)), fBranches(bl), fPrevData(pd), fValues(fNSlots)
647  {
648  }
649 
650  RFilter(const RFilter &) = delete;
651  RFilter &operator=(const RFilter &) = delete;
652 
653  bool CheckFilters(unsigned int slot, Long64_t entry) final
654  {
655  if (entry != fLastCheckedEntry[slot]) {
656  if (!fPrevData.CheckFilters(slot, entry)) {
657  // a filter upstream returned false, cache the result
658  fLastResult[slot] = false;
659  } else {
660  // evaluate this filter, cache the result
661  auto passed = CheckFilterHelper(slot, entry, TypeInd_t());
662  passed ? ++fAccepted[slot] : ++fRejected[slot];
663  fLastResult[slot] = passed;
664  }
665  fLastCheckedEntry[slot] = entry;
666  }
667  return fLastResult[slot];
668  }
669 
670  template <std::size_t... S>
671  bool CheckFilterHelper(unsigned int slot, Long64_t entry, std::index_sequence<S...>)
672  {
673  return fFilter(std::get<S>(fValues[slot]).Get(entry)...);
674  // silence "unused parameter" warnings in gcc
675  (void)slot;
676  (void)entry;
677  }
678 
679  void InitSlot(TTreeReader *r, unsigned int slot) final
680  {
681  RDFInternal::InitRDFValues(slot, fValues[slot], r, fBranches, fLoopManager->GetCustomColumnNames(),
682  fLoopManager->GetBookedColumns(), TypeInd_t());
683  }
684 
685  // recursive chain of `Report`s
686  void Report(ROOT::RDF::RCutFlowReport &rep) const final { PartialReport(rep); }
687 
689  {
690  fPrevData.PartialReport(rep);
691  FillReport(rep);
692  }
693 
694  void StopProcessing() final
695  {
696  ++fNStopsReceived;
697  if (fNStopsReceived == fNChildren)
698  fPrevData.StopProcessing();
699  }
700 
701  void IncrChildrenCount() final
702  {
703  ++fNChildren;
704  // propagate "children activation" upstream. named filters do the propagation via `TriggerChildrenCount`.
705  if (fNChildren == 1 && fName.empty())
706  fPrevData.IncrChildrenCount();
707  }
708 
709  void TriggerChildrenCount() final
710  {
711  assert(!fName.empty()); // this method is to only be called on named filters
712  fPrevData.IncrChildrenCount();
713  }
714 
715  virtual void ClearValueReaders(unsigned int slot) final
716  {
717  RDFInternal::ResetRDFValueTuple(fValues[slot], TypeInd_t());
718  }
719 };
720 
721 class RRangeBase {
722 protected:
723  RLoopManager *fLoopManager; ///< A raw pointer to the RLoopManager at the root of this functional graph. It is only
724  /// guaranteed to contain a valid address during an event loop.
725  unsigned int fStart;
726  unsigned int fStop;
727  unsigned int fStride;
728  Long64_t fLastCheckedEntry{-1};
729  bool fLastResult{true};
730  ULong64_t fNProcessedEntries{0};
731  unsigned int fNChildren{0}; ///< Number of nodes of the functional graph hanging from this object
732  unsigned int fNStopsReceived{0}; ///< Number of times that a children node signaled to stop processing entries.
733  bool fHasStopped{false}; ///< True if the end of the range has been reached
734  const unsigned int fNSlots; ///< Number of thread slots used by this node, inherited from parent node.
735 
736  void ResetCounters();
737 
738 public:
739  RRangeBase(RLoopManager *implPtr, unsigned int start, unsigned int stop, unsigned int stride,
740  const unsigned int nSlots);
741  RRangeBase &operator=(const RRangeBase &) = delete;
742  virtual ~RRangeBase() = default;
743 
744  RLoopManager *GetLoopManagerUnchecked() const;
745  virtual bool CheckFilters(unsigned int slot, Long64_t entry) = 0;
746  virtual void Report(ROOT::RDF::RCutFlowReport &) const = 0;
747  virtual void PartialReport(ROOT::RDF::RCutFlowReport &) const = 0;
748  virtual void IncrChildrenCount() = 0;
749  virtual void StopProcessing() = 0;
751  {
752  fNChildren = 0;
753  fNStopsReceived = 0;
754  }
755  void InitNode() { ResetCounters(); }
756 };
757 
758 template <typename PrevData>
759 class RRange final : public RRangeBase {
760  PrevData &fPrevData;
761 
762 public:
763  RRange(unsigned int start, unsigned int stop, unsigned int stride, PrevData &pd)
764  : RRangeBase(pd.GetLoopManagerUnchecked(), start, stop, stride, pd.GetLoopManagerUnchecked()->GetNSlots()),
765  fPrevData(pd)
766  {
767  }
768 
769  RRange(const RRange &) = delete;
770  RRange &operator=(const RRange &) = delete;
771 
772  /// Ranges act as filters when it comes to selecting entries that downstream nodes should process
773  bool CheckFilters(unsigned int slot, Long64_t entry) final
774  {
775  if (entry != fLastCheckedEntry) {
776  if (fHasStopped)
777  return false;
778  if (!fPrevData.CheckFilters(slot, entry)) {
779  // a filter upstream returned false, cache the result
780  fLastResult = false;
781  } else {
782  // apply range filter logic, cache the result
783  ++fNProcessedEntries;
784  if (fNProcessedEntries <= fStart || (fStop > 0 && fNProcessedEntries > fStop) ||
785  (fStride != 1 && fNProcessedEntries % fStride != 0))
786  fLastResult = false;
787  else
788  fLastResult = true;
789  if (fNProcessedEntries == fStop) {
790  fHasStopped = true;
791  fPrevData.StopProcessing();
792  }
793  }
794  fLastCheckedEntry = entry;
795  }
796  return fLastResult;
797  }
798 
799  // recursive chain of `Report`s
800  // RRange simply forwards these calls to the previous node
801  void Report(ROOT::RDF::RCutFlowReport &rep) const final { fPrevData.PartialReport(rep); }
802 
803  void PartialReport(ROOT::RDF::RCutFlowReport &rep) const final { fPrevData.PartialReport(rep); }
804 
805  void StopProcessing() final
806  {
807  ++fNStopsReceived;
808  if (fNStopsReceived == fNChildren && !fHasStopped)
809  fPrevData.StopProcessing();
810  }
811 
812  void IncrChildrenCount() final
813  {
814  ++fNChildren;
815  // propagate "children activation" upstream
816  if (fNChildren == 1)
817  fPrevData.IncrChildrenCount();
818  }
819 };
820 
821 } // namespace RDF
822 } // namespace Detail
823 
824 // method implementations
825 namespace Internal {
826 namespace RDF {
827 
828 template <typename T, bool B>
830 {
831  fCustomColumns.emplace_back(customColumn);
832  if (customColumn->GetTypeId() != typeid(T))
833  throw std::runtime_error(
834  std::string("TColumnValue: type specified for column \"" + customColumn->GetName() + "\" is ") +
835  TypeID2TypeName(typeid(T)) + " but temporary column has type " + TypeID2TypeName(customColumn->GetTypeId()));
836 
837  if (customColumn->IsDataSourceColumn()) {
838  fColumnKind = EColumnKind::kDataSource;
839  fDSValuePtrs.emplace_back(static_cast<T **>(customColumn->GetValuePtr(slot)));
840  } else {
841  fColumnKind = EColumnKind::kCustomColumn;
842  fCustomValuePtrs.emplace_back(static_cast<T *>(customColumn->GetValuePtr(slot)));
843  }
844  fSlot = slot;
845 }
846 
847 // This method is executed inside the event-loop, many times per entry
848 // If need be, the if statement can be avoided using thunks
849 // (have both branches inside functions and have a pointer to the branch to be executed)
850 template <typename T, bool B>
851 template <typename U, typename std::enable_if<!TColumnValue<U>::fgMustUseRVec, int>::type>
853 {
854  if (fColumnKind == EColumnKind::kTree) {
855  return *(fTreeReaders.back()->Get());
856  } else {
857  fCustomColumns.back()->Update(fSlot, entry);
858  return fColumnKind == EColumnKind::kCustomColumn ? *fCustomValuePtrs.back() : **fDSValuePtrs.back();
859  }
860 }
861 
862 /// This overload is used to return arrays (i.e. types that are read into a RVec)
863 template <typename T, bool B>
864 template <typename U, typename std::enable_if<TColumnValue<U>::fgMustUseRVec, int>::type>
866 {
867  if (fColumnKind == EColumnKind::kTree) {
868  auto &readerArray = *fTreeReaders.back();
869  // We only use TTreeReaderArrays to read columns that users flagged as type `RVec`, so we need to check
870  // that the branch stores the array as contiguous memory that we can actually wrap in an `RVec`.
871  // Currently we need the first entry to have been loaded to perform the check
872  // TODO Move check to `MakeProxy` once Axel implements this kind of check in TTreeReaderArray using
873  // TBranchProxy
874 
875  if (EStorageType::kUnknown == fStorageType && readerArray.GetSize() > 1) {
876  // We can decide since the array is long enough
877  fStorageType = (1 == (&readerArray[1] - &readerArray[0])) ? EStorageType::kContiguous : EStorageType::kSparse;
878  }
879 
880  const auto readerArraySize = readerArray.GetSize();
881  if (EStorageType::kContiguous == fStorageType ||
882  (EStorageType::kUnknown == fStorageType && readerArray.GetSize() < 2)) {
883  if (readerArraySize > 0) {
884  // trigger loading of the contens of the TTreeReaderArray
885  // the address of the first element in the reader array is not necessarily equal to
886  // the address returned by the GetAddress method
887  auto readerArrayAddr = &readerArray.At(0);
888  T tvec(readerArrayAddr, readerArraySize);
889  swap(fRVec, tvec);
890  } else {
891  T emptyVec{};
892  swap(fRVec, emptyVec);
893  }
894  } else {
895  // The storage is not contiguous or we don't know yet: we cannot but copy into the tvec
896 #ifndef NDEBUG
897  if (!fCopyWarningPrinted) {
898  Warning("TColumnValue::Get", "Branch %s hangs from a non-split branch. A copy is being performed in order "
899  "to properly read the content.",
900  readerArray.GetBranchName());
901  fCopyWarningPrinted = true;
902  }
903 #else
904  (void)fCopyWarningPrinted;
905 #endif
906  if (readerArraySize > 0) {
907  (void)readerArray.At(0); // trigger deserialisation
908  T tvec(readerArray.begin(), readerArray.end());
909  swap(fRVec, tvec);
910  } else {
911  T emptyVec{};
912  swap(fRVec, emptyVec);
913  }
914  }
915  return fRVec;
916 
917  } else {
918  fCustomColumns.back()->Update(fSlot, entry);
919  return fColumnKind == EColumnKind::kCustomColumn ? *fCustomValuePtrs.back() : **fDSValuePtrs.back();
920  }
921 }
922 
923 } // namespace RDF
924 } // namespace Internal
925 } // namespace ROOT
926 #endif // ROOT_RDFNODES
typename std::conditional< MustUseRVec, TTreeReaderArray< ColumnValue_t >, TTreeReaderValue< ColumnValue_t > >::type TreeReader_t
Definition: RDFNodes.hxx:247
std::vector< RDFInternal::RDFValueTuple_t< ColumnTypes_t > > fValues
Definition: RDFNodes.hxx:641
std::string GetName(const std::string &scope_name)
Definition: Cppyy.cxx:145
const unsigned int fNSlots
Number of thread slots used by this node.
Definition: RDFNodes.hxx:343
TSlotStack(unsigned int size)
Definition: RDFNodes.hxx:65
void AddCustomColumnName(std::string_view name)
Definition: RDFNodes.hxx:209
long long Long64_t
Definition: RtypesCore.h:69
void Run(unsigned int slot, Long64_t entry) final
Definition: RDFNodes.hxx:390
TTreeReader is a simple, robust and fast interface to read values from a TTree, TChain or TNtuple...
Definition: TTreeReader.h:43
Namespace for new ROOT classes and functions.
Definition: StringConv.hxx:21
std::shared_ptr< RFilterBase > FilterBasePtr_t
Definition: RDFNodes.hxx:83
virtual void * GetValuePtr(unsigned int slot)=0
std::shared_ptr< RDFInternal::RActionBase > ActionBasePtr_t
Definition: RDFNodes.hxx:78
RFilter(FilterF &&f, const ColumnNames_t &bl, PrevDataFrame &pd, std::string_view name="")
Definition: RDFNodes.hxx:644
std::vector< FilterBasePtr_t > FilterBaseVec_t
Definition: RDFNodes.hxx:84
typename TRDFValueTuple< BranchType >::type RDFValueTuple_t
Definition: RDFNodes.hxx:327
PrevDataFrame & fPrevData
Definition: RDFNodes.hxx:640
void swap(TDirectoryEntry &e1, TDirectoryEntry &e2) noexcept
bool CheckFilters(unsigned int slot, Long64_t entry) final
Ranges act as filters when it comes to selecting entries that downstream nodes should process...
Definition: RDFNodes.hxx:773
double T(double x)
Definition: ChebyshevPol.h:34
std::map< std::string, std::string > fAliasColumnNameMap
ColumnNameAlias-columnName pairs.
Definition: RDFNodes.hxx:153
std::vector< ActionBasePtr_t > ActionBaseVec_t
Definition: RDFNodes.hxx:79
void PartialReport(ROOT::RDF::RCutFlowReport &rep) const final
Definition: RDFNodes.hxx:803
RLoopManager * fLoopManager
A raw pointer to the RLoopManager at the root of this functional graph.
Definition: RDFNodes.hxx:564
fill
Definition: fit1_py.py:6
FilterBaseVec_t fBookedNamedFilters
Contains a subset of fBookedFilters, i.e. only the named filters.
Definition: RDFNodes.hxx:134
ActionBaseVec_t fBookedActions
Definition: RDFNodes.hxx:132
void UpdateHelper(unsigned int slot, Long64_t entry, std::index_sequence< S... >, TypeList< BranchTypes... >, TCCHelperTypes::TSlotAndEntry *)
Definition: RDFNodes.hxx:550
void SetTmpColumn(unsigned int slot, RCustomColumnBase *tmpColumn)
Definition: RDFNodes.hxx:829
void UpdateHelper(unsigned int slot, Long64_t entry, std::index_sequence< S... >, TypeList< BranchTypes... >, TCCHelperTypes::TSlot *)
Definition: RDFNodes.hxx:540
#define f(i)
Definition: RSha256.hxx:104
const ColumnNames_t fDefaultColumns
Definition: RDFNodes.hxx:143
std::vector< T * > fCustomValuePtrs
Non-owning ptrs to the value of a custom column.
Definition: RDFNodes.hxx:264
typename RDFInternal::RemoveFirstTwoParametersIf< std::is_same< TSlotAndEntry, UHT_t >::value, BranchTypesTmp_t >::type ColumnTypes_t
Definition: RDFNodes.hxx:483
unsigned int GetNSlots()
Definition: RDFUtils.cxx:244
STL namespace.
std::make_index_sequence< ColumnTypes_t::list_size > TypeInd_t
Definition: RDFNodes.hxx:363
void ClearValueReaders(unsigned int slot)
Definition: RDFNodes.hxx:412
void * PartialUpdate(unsigned int slot) final
This method is invoked to update a partial result during the event loop, right before passing the res...
Definition: RDFNodes.hxx:417
A spin mutex class which respects the STL interface for mutexes.
Definition: TSpinMutex.hxx:40
std::string TypeID2TypeName(const std::type_info &id)
Returns the name of a type starting from its type_info An empty string is returned in case of failure...
Definition: RDFUtils.cxx:82
TCallback(ULong64_t everyN, Callback_t &&f, unsigned int nSlots)
Definition: RDFNodes.hxx:100
const ELoopType fLoopType
The kind of event loop that is going to be run (e.g. on ROOT files, on no files)
Definition: RDFNodes.hxx:149
typename CallableTraits< FilterF >::arg_types ColumnTypes_t
Definition: RDFNodes.hxx:635
std::map< std::string, RCustomColumnBasePtr_t > fBookedCustomColumns
Definition: RDFNodes.hxx:135
RLoopManager * fLoopManager
A raw pointer to the RLoopManager at the root of this functional graph.
Definition: RDFNodes.hxx:340
RAction(Helper &&h, const ColumnNames_t &bl, PrevDataFrame &pd)
Definition: RDFNodes.hxx:371
std::vector< RDFValueTuple_t< ColumnTypes_t > > fValues
Definition: RDFNodes.hxx:368
typename CallableTraits< F >::ret_type ret_type
Definition: RDFNodes.hxx:485
void IncrChildrenCount() final
Definition: RDFNodes.hxx:812
void * GetValuePtr(unsigned int slot) final
Definition: RDFNodes.hxx:513
void MakeProxy(TTreeReader *r, const std::string &bn)
Definition: RDFNodes.hxx:285
void AddDataSourceColumn(std::string_view name)
Definition: RDFNodes.hxx:207
void ReturnSlot(unsigned int slotNumber)
Definition: RDFNodes.cxx:191
const ColumnNames_t fBranches
Definition: RDFNodes.hxx:639
unsigned int GetNSlots() const
Definition: RDFNodes.hxx:197
ULong64_t GetNEmptyEntries() const
Definition: RDFNodes.hxx:189
void TriggerChildrenCount() final
Definition: RDFNodes.hxx:709
RLoopManager * fLoopManager
A raw pointer to the RLoopManager at the root of this functional graph.
Definition: RDFNodes.hxx:723
std::vector< RDFInternal::RDFValueTuple_t< ColumnTypes_t > > fValues
Definition: RDFNodes.hxx:494
void function(const Char_t *name_, T fun, const Char_t *docstring=0)
Definition: RExports.h:146
A wrapper around a concrete RFilter, which forwards all calls to it RJittedFilter is the type of the ...
Definition: RDFNodes.hxx:611
std::vector< unsigned int > fBuf
Definition: RDFNodes.hxx:60
const ColumnNames_t fBranches
Definition: RDFNodes.hxx:491
T & Get(Long64_t entry)
This overload is used to return scalar quantities (i.e. types that are not read into a RVec) ...
Definition: RDFNodes.hxx:852
virtual void ResetReportCount()
Definition: RDFNodes.hxx:595
Extracts data from a TTree.
EColumnKind
TColumnValue has a slightly different behaviour whether the column comes from a TTreeReader, a RDataFrame Define or a RDataSource.
Definition: RDFNodes.hxx:251
std::vector< ULong64_t > fCounters
Definition: RDFNodes.hxx:97
unsigned int GetID() const
Definition: RDFNodes.hxx:212
void FinalizeSlot(unsigned int slot) final
Definition: RDFNodes.hxx:406
TOneTimeCallback(Callback_t &&f, unsigned int nSlots)
Definition: RDFNodes.hxx:121
TCCHelperTypes::TSlot TSlot
Definition: RDFNodes.hxx:475
RooArgSet S(const RooAbsArg &v1)
#define F(x, y, z)
const std::type_info & GetTypeId() const
Definition: RDFNodes.hxx:524
const ColumnNames_t & GetCustomColumnNames() const
Definition: RDFNodes.hxx:185
RJittedFilter(RLoopManager *lm, std::string_view name)
Definition: RDFNodes.hxx:615
ROOT::R::TRInterface & r
Definition: Object.C:4
void Initialize(Bool_t useTMVAStyle=kTRUE)
Definition: tmvaglob.cxx:176
Helper class that updates and returns TTree branches as well as RDataFrame temporary columns...
Definition: RDFNodes.hxx:243
std::vector< TCallback > fCallbacks
Registered callbacks.
Definition: RDFNodes.hxx:154
std::shared_ptr< RCustomColumnBase > RCustomColumnBasePtr_t
Definition: RDFNodes.hxx:81
virtual const std::type_info & GetTypeId() const =0
RCustomColumn(std::string_view name, F &&expression, const ColumnNames_t &bl, RLoopManager *lm, bool isDSColumn=false)
Definition: RDFNodes.hxx:497
bool CheckFilters(unsigned int slot, Long64_t entry) final
Definition: RDFNodes.hxx:653
void Update(unsigned int slot, Long64_t entry) final
Definition: RDFNodes.hxx:515
auto PartialUpdateImpl(unsigned int slot) -> decltype(std::declval< H >().PartialUpdate(slot),(void *)(nullptr))
Definition: RDFNodes.hxx:423
void Report(ROOT::RDF::RCutFlowReport &rep) const final
Definition: RDFNodes.hxx:801
std::make_index_sequence< ColumnTypes_t::list_size > TypeInd_t
Definition: RDFNodes.hxx:484
std::vector< std::shared_ptr< bool > > fResProxyReadiness
Definition: RDFNodes.hxx:138
void UpdateHelper(unsigned int slot, Long64_t entry, std::index_sequence< S... >, TypeList< BranchTypes... >, TCCHelperTypes::TNothing *)
Definition: RDFNodes.hxx:530
std::string fToJit
string containing all BuildAndBook actions that should be jitted before running
Definition: RDFNodes.hxx:150
std::vector< Long64_t > fLastCheckedEntry
Definition: RDFNodes.hxx:446
void InitSlot(TTreeReader *r, unsigned int slot) final
Definition: RDFNodes.hxx:383
virtual void ResetChildrenCount()
Definition: RDFNodes.hxx:589
void StopProcessing() final
Definition: RDFNodes.hxx:694
std::tuple< TColumnValue< BranchTypes >... > type
Definition: RDFNodes.hxx:323
void Reset(Detail::TBranchProxy *x)
void Warning(const char *location, const char *msgfmt,...)
void AddColumnAlias(const std::string &alias, const std::string &colName)
Definition: RDFNodes.hxx:208
#define h(i)
Definition: RSha256.hxx:106
Lightweight storage for a collection of types.
Definition: TypeTraits.hxx:27
std::make_index_sequence< ColumnTypes_t::list_size > TypeInd_t
Definition: RDFNodes.hxx:636
typename CallableTraits< F >::arg_types FunParamTypes_t
Definition: RDFNodes.hxx:479
bool CheckFilterHelper(unsigned int slot, Long64_t entry, std::index_sequence< S... >)
Definition: RDFNodes.hxx:671
PrevDataFrame & fPrevData
Definition: RDFNodes.hxx:367
const std::map< std::string, RCustomColumnBasePtr_t > & GetBookedColumns() const
Definition: RDFNodes.hxx:187
std::vector< RCustomColumnBase * > fCustomColumns
Non-owning ptrs to the node responsible for the custom column. Needed when querying custom values...
Definition: RDFNodes.hxx:268
Describe directory structure in memory.
Definition: TDirectory.h:34
int type
Definition: TGX11.cxx:120
unsigned long long ULong64_t
Definition: RtypesCore.h:70
basic_string_view< char > string_view
Definition: RStringView.hxx:35
ROOT type_traits extensions.
Definition: TypeTraits.hxx:23
std::shared_ptr< RRangeBase > RangeBasePtr_t
Definition: RDFNodes.hxx:86
static constexpr double s
void Exec(unsigned int slot, Long64_t entry, std::index_sequence< S... >)
Definition: RDFNodes.hxx:398
virtual void ClearValueReaders(unsigned int slot) final
Definition: RDFNodes.hxx:715
const ColumnNames_t fBranches
Definition: RDFNodes.hxx:366
void Report(ROOT::RDF::RCutFlowReport &rep) const final
Definition: RDFNodes.hxx:686
void ClearValueReaders(unsigned int slot) final
Definition: RDFNodes.hxx:559
Binding & operator=(OUT(*fun)(void))
ColumnNames_t fDefinedDataSourceColumns
List of data-source columns that have been Defined so far.
Definition: RDFNodes.hxx:152
const std::unique_ptr< RDataSource > fDataSource
Owning pointer to a data-source object. Null if no data-source.
Definition: RDFNodes.hxx:151
typename RDFInternal::RemoveFirstParameterIf< std::is_same< TSlot, UHT_t >::value, FunParamTypes_t >::type BranchTypesTmp_t
Definition: RDFNodes.hxx:481
typename std::conditional< std::is_same< ret_type, bool >::value, std::deque< ret_type >, std::vector< ret_type > >::type ValuesPerSlot_t
Definition: RDFNodes.hxx:488
typedef void((*Func_t)())
std::vector< Long64_t > fLastCheckedEntry
Definition: RDFNodes.hxx:566
void SetTree(const std::shared_ptr< TTree > &tree)
Definition: RDFNodes.hxx:202
void ToJit(const std::string &s)
Definition: RDFNodes.hxx:205
RDataSource * GetDataSource() const
Definition: RDFNodes.hxx:190
void InitRDFValues(unsigned int slot, RDFValueTuple &valueTuple, TTreeReader *r, const ColumnNames_t &bn, const ColumnNames_t &tmpbn, const std::map< std::string, std::shared_ptr< RCustomColumnBase >> &customCols, std::index_sequence< S... >)
Initialize a tuple of TColumnValues.
void IncrChildrenCount() final
Definition: RDFNodes.hxx:701
const unsigned int fNSlots
number of thread slots used by this node, inherited from parent node.
Definition: RDFNodes.hxx:444
EStorageType
Enumerator for the different properties of the branch storage in memory.
Definition: RDFNodes.hxx:270
void PartialReport(ROOT::RDF::RCutFlowReport &) const
End of recursive chain of calls, does nothing.
Definition: RDFNodes.hxx:201
RRange(unsigned int start, unsigned int stop, unsigned int stride, PrevData &pd)
Definition: RDFNodes.hxx:763
#define c(i)
Definition: RSha256.hxx:101
const unsigned int fNSlots
Number of thread slots used by this node, inherited from parent node.
Definition: RDFNodes.hxx:573
Definition: tree.py:1
const bool fIsDataSourceColumn
does the custom column refer to a data-source column? (or a user-define column?)
Definition: RDFNodes.hxx:445
RVec< ColumnValue_t > fRVec
If MustUseRVec, i.e. we are reading an array, we return a reference to this RVec to clients...
Definition: RDFNodes.hxx:275
const std::map< std::string, std::string > & GetAliasMap() const
Definition: RDFNodes.hxx:210
RLoopManager * fLoopManager
A raw pointer to the RLoopManager at the root of this functional graph.
Definition: RDFNodes.hxx:439
FilterBaseVec_t fBookedFilters
Definition: RDFNodes.hxx:133
void InitSlot(TTreeReader *r, unsigned int slot) final
Definition: RDFNodes.hxx:507
void StopProcessing() final
Definition: RDFNodes.hxx:805
RDataSource defines an API that RDataFrame can use to read arbitrary data formats.
Definition: RDataSource.hxx:91
void ResetRDFValueTuple(ValueTuple &values, std::index_sequence< S... >)
Clear the proxies of a tuple of TColumnValues.
Definition: RDFNodes.hxx:331
std::function< void(unsigned int)> Callback_t
Definition: RDFNodes.hxx:93
const unsigned int fNSlots
Number of thread slots used by this node, inherited from parent node.
Definition: RDFNodes.hxx:734
const ColumnNames_t & GetDefinedDataSourceColumns() const
Definition: RDFNodes.hxx:206
TCCHelperTypes::TSlotAndEntry TSlotAndEntry
Definition: RDFNodes.hxx:476
void InitSlot(TTreeReader *r, unsigned int slot) final
Definition: RDFNodes.hxx:679
make_integer_sequence< size_t, _Np > make_index_sequence
std::vector< std::unique_ptr< TreeReader_t > > fTreeReaders
Owning ptrs to a TTreeReaderValue or TTreeReaderArray. Only used for Tree columns.
Definition: RDFNodes.hxx:262
char name[80]
Definition: TGX11.cxx:109
void TriggerChildrenCount() final
Definition: RDFNodes.hxx:404
void PartialReport(ROOT::RDF::RCutFlowReport &rep) const final
Definition: RDFNodes.hxx:688
std::vector< T ** > fDSValuePtrs
Non-owning ptrs to the value of a data-source column.
Definition: RDFNodes.hxx:266
TCCHelperTypes::TNothing TNothing
Definition: RDFNodes.hxx:474
std::vector< RangeBasePtr_t > RangeBaseVec_t
Definition: RDFNodes.hxx:87
typename std::conditional< MustUseRVec, TakeFirstParameter_t< T >, T >::type ColumnValue_t
Definition: RDFNodes.hxx:245
std::vector< TOneTimeCallback > fCallbacksOnce
Registered callbacks to invoke just once before running the loop.
Definition: RDFNodes.hxx:155
ColumnNames_t fCustomColumnNames
Contains names of all custom columns defined in the functional graph.
Definition: RDFNodes.hxx:136