Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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
17#include <ROOT/RDF/RDefine.hxx>
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#include <TROOT.h> // IsImplicitMTEnabled
31
32#include <deque>
33#include <functional>
34#include <map>
35#include <memory>
36#include <string>
37#include <type_traits>
38#include <typeinfo>
39#include <vector>
40#include <unordered_map>
41
42class TObjArray;
43class TTree;
44namespace ROOT {
45namespace Detail {
46namespace RDF {
47class RNodeBase;
48}
49}
50namespace RDF {
51template <typename T>
52class RResultPtr;
53template<typename T, typename V>
54class RInterface;
56class RDataSource;
57} // namespace RDF
58
59} // namespace ROOT
60
61/// \cond HIDDEN_SYMBOLS
62
63namespace ROOT {
64namespace Internal {
65namespace RDF {
66using namespace ROOT::Detail::RDF;
67using namespace ROOT::RDF;
68namespace TTraits = ROOT::TypeTraits;
70
72
73std::string DemangleTypeIdName(const std::type_info &typeInfo);
74
76ConvertRegexToColumns(const ColumnNames_t &colNames, std::string_view columnNameRegexp, std::string_view callerName);
77
78/// An helper object that sets and resets gErrorIgnoreLevel via RAII.
79class RIgnoreErrorLevelRAII {
80private:
81 int fCurIgnoreErrorLevel = gErrorIgnoreLevel;
82
83public:
84 RIgnoreErrorLevelRAII(int errorIgnoreLevel) { gErrorIgnoreLevel = errorIgnoreLevel; }
85 ~RIgnoreErrorLevelRAII() { gErrorIgnoreLevel = fCurIgnoreErrorLevel; }
86};
87
88/****** BuildAction overloads *******/
89
90// clang-format off
91/// This namespace defines types to be used for tag dispatching in RInterface.
92namespace ActionTags {
93struct Histo1D{};
94struct Histo2D{};
95struct Histo3D{};
96struct Graph{};
97struct Profile1D{};
98struct Profile2D{};
99struct Min{};
100struct Max{};
101struct Sum{};
102struct Mean{};
103struct Fill{};
104struct StdDev{};
105struct Display{};
106struct Snapshot{};
107}
108// clang-format on
109
110template <typename T, bool ISV6HISTO = std::is_base_of<TH1, T>::value>
111struct HistoUtils {
112 static void SetCanExtendAllAxes(T &h) { h.SetCanExtend(::TH1::kAllAxes); }
113 static bool HasAxisLimits(T &h)
114 {
115 auto xaxis = h.GetXaxis();
116 return !(xaxis->GetXmin() == 0. && xaxis->GetXmax() == 0.);
117 }
118};
119
120template <typename T>
121struct HistoUtils<T, false> {
122 static void SetCanExtendAllAxes(T &) {}
123 static bool HasAxisLimits(T &) { return true; }
124};
125
126// Generic filling (covers Histo2D, Histo3D, Profile1D and Profile2D actions, with and without weights)
127template <typename... ColTypes, typename ActionTag, typename ActionResultType, typename PrevNodeType>
128std::unique_ptr<RActionBase>
129BuildAction(const ColumnNames_t &bl, const std::shared_ptr<ActionResultType> &h, const unsigned int nSlots,
130 std::shared_ptr<PrevNodeType> prevNode, ActionTag, const RDFInternal::RBookedDefines &defines)
131{
132 using Helper_t = FillParHelper<ActionResultType>;
133 using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<ColTypes...>>;
134 return std::make_unique<Action_t>(Helper_t(h, nSlots), bl, std::move(prevNode), defines);
135}
136
137// Histo1D filling (must handle the special case of distinguishing FillParHelper and FillHelper
138template <typename... ColTypes, typename PrevNodeType>
139std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<::TH1D> &h,
140 const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
141 ActionTags::Histo1D, const RDFInternal::RBookedDefines &defines)
142{
143 auto hasAxisLimits = HistoUtils<::TH1D>::HasAxisLimits(*h);
144
145 if (hasAxisLimits) {
146 using Helper_t = FillParHelper<::TH1D>;
147 using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<ColTypes...>>;
148 return std::make_unique<Action_t>(Helper_t(h, nSlots), bl, std::move(prevNode), defines);
149 } else {
150 using Helper_t = FillHelper;
151 using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<ColTypes...>>;
152 return std::make_unique<Action_t>(Helper_t(h, nSlots), bl, std::move(prevNode), defines);
153 }
154}
155
156template <typename... ColTypes, typename PrevNodeType>
157std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<TGraph> &g,
158 const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
159 ActionTags::Graph, const RDFInternal::RBookedDefines &defines)
160{
161 using Helper_t = FillTGraphHelper;
162 using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<ColTypes...>>;
163 return std::make_unique<Action_t>(Helper_t(g, nSlots), bl, std::move(prevNode), defines);
164}
165
166// Min action
167template <typename ColType, typename PrevNodeType, typename ActionResultType>
168std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<ActionResultType> &minV,
169 const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
170 ActionTags::Min, const RDFInternal::RBookedDefines &defines)
171{
172 using Helper_t = MinHelper<ActionResultType>;
173 using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<ColType>>;
174 return std::make_unique<Action_t>(Helper_t(minV, nSlots), bl, std::move(prevNode), defines);
175}
176
177// Max action
178template <typename ColType, typename PrevNodeType, typename ActionResultType>
179std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<ActionResultType> &maxV,
180 const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
181 ActionTags::Max, const RDFInternal::RBookedDefines &defines)
182{
183 using Helper_t = MaxHelper<ActionResultType>;
184 using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<ColType>>;
185 return std::make_unique<Action_t>(Helper_t(maxV, nSlots), bl, std::move(prevNode), defines);
186}
187
188// Sum action
189template <typename ColType, typename PrevNodeType, typename ActionResultType>
190std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<ActionResultType> &sumV,
191 const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
192 ActionTags::Sum, const RDFInternal::RBookedDefines &defines)
193{
194 using Helper_t = SumHelper<ActionResultType>;
195 using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<ColType>>;
196 return std::make_unique<Action_t>(Helper_t(sumV, nSlots), bl, std::move(prevNode), defines);
197}
198
199// Mean action
200template <typename ColType, typename PrevNodeType>
201std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<double> &meanV,
202 const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
203 ActionTags::Mean, const RDFInternal::RBookedDefines &defines)
204{
205 using Helper_t = MeanHelper;
206 using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<ColType>>;
207 return std::make_unique<Action_t>(Helper_t(meanV, nSlots), bl, std::move(prevNode), defines);
208}
209
210// Standard Deviation action
211template <typename ColType, typename PrevNodeType>
212std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<double> &stdDeviationV,
213 const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode,
214 ActionTags::StdDev, const RDFInternal::RBookedDefines &defines)
215{
216 using Helper_t = StdDevHelper;
217 using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<ColType>>;
218 return std::make_unique<Action_t>(Helper_t(stdDeviationV, nSlots), bl, prevNode, defines);
219}
220
221// Display action
222template <typename... ColTypes, typename PrevNodeType>
223std::unique_ptr<RActionBase> BuildAction(const ColumnNames_t &bl, const std::shared_ptr<RDisplay> &d,
224 const unsigned int, std::shared_ptr<PrevNodeType> prevNode,
225 ActionTags::Display, const RDFInternal::RBookedDefines &defines)
226{
227 using Helper_t = DisplayHelper<PrevNodeType>;
228 using Action_t = RAction<Helper_t, PrevNodeType, TTraits::TypeList<ColTypes...>>;
229 return std::make_unique<Action_t>(Helper_t(d, prevNode), bl, prevNode, defines);
230}
231
232struct SnapshotHelperArgs {
233 std::string fFileName;
234 std::string fDirName;
235 std::string fTreeName;
236 std::vector<std::string> fOutputColNames;
238};
239
240// Snapshot action
241template <typename... ColTypes, typename PrevNodeType>
242std::unique_ptr<RActionBase>
243BuildAction(const ColumnNames_t &colNames, const std::shared_ptr<SnapshotHelperArgs> &snapHelperArgs,
244 const unsigned int nSlots, std::shared_ptr<PrevNodeType> prevNode, ActionTags::Snapshot,
245 const RDFInternal::RBookedDefines &defines)
246{
247 const auto &filename = snapHelperArgs->fFileName;
248 const auto &dirname = snapHelperArgs->fDirName;
249 const auto &treename = snapHelperArgs->fTreeName;
250 const auto &outputColNames = snapHelperArgs->fOutputColNames;
251 const auto &options = snapHelperArgs->fOptions;
252
253 std::unique_ptr<RActionBase> actionPtr;
255 // single-thread snapshot
256 using Helper_t = SnapshotHelper<ColTypes...>;
257 using Action_t = RAction<Helper_t, PrevNodeType>;
258 actionPtr.reset(new Action_t(Helper_t(filename, dirname, treename, colNames, outputColNames, options), colNames,
259 prevNode, defines));
260 } else {
261 // multi-thread snapshot
262 using Helper_t = SnapshotHelperMT<ColTypes...>;
263 using Action_t = RAction<Helper_t, PrevNodeType>;
264 actionPtr.reset(new Action_t(Helper_t(nSlots, filename, dirname, treename, colNames, outputColNames, options),
265 colNames, prevNode, defines));
266 }
267 return actionPtr;
268}
269
270/****** end BuildAndBook ******/
271
272template <typename Filter>
273void CheckFilter(Filter &)
274{
275 using FilterRet_t = typename RDF::CallableTraits<Filter>::ret_type;
276 static_assert(std::is_convertible<FilterRet_t, bool>::value,
277 "filter expression returns a type that is not convertible to bool");
278}
279
280void CheckDefine(std::string_view definedCol, TTree *treePtr, const ColumnNames_t &customCols,
281 const std::map<std::string, std::string> &aliasMap, const ColumnNames_t &dataSourceColumns);
282
283std::string PrettyPrintAddr(const void *const addr);
284
285void BookFilterJit(const std::shared_ptr<RJittedFilter> &jittedFilter, std::shared_ptr<RNodeBase> *prevNodeOnHeap,
286 std::string_view name, std::string_view expression,
287 const std::map<std::string, std::string> &aliasMap, const ColumnNames_t &branches,
288 const RDFInternal::RBookedDefines &customCols, TTree *tree, RDataSource *ds);
289
290std::shared_ptr<RJittedDefine> BookDefineJit(std::string_view name, std::string_view expression, RLoopManager &lm,
291 RDataSource *ds, const RDFInternal::RBookedDefines &customCols,
292 const ColumnNames_t &branches,
293 std::shared_ptr<RNodeBase> *prevNodeOnHeap);
294
295std::string JitBuildAction(const ColumnNames_t &bl, std::shared_ptr<RDFDetail::RNodeBase> *prevNode,
296 const std::type_info &art, const std::type_info &at, void *rOnHeap, TTree *tree,
297 const unsigned int nSlots, const RDFInternal::RBookedDefines &defines,
298 RDataSource *ds, std::weak_ptr<RJittedAction> *jittedActionOnHeap);
299
300// Allocate a weak_ptr on the heap, return a pointer to it. The user is responsible for deleting this weak_ptr.
301// This function is meant to be used by RInterface's methods that book code for jitting.
302// The problem it solves is that we generate code to be lazily jitted with the addresses of certain objects in them,
303// and we need to check those objects are still alive when the generated code is finally jitted and executed.
304// So we pass addresses to weak_ptrs allocated on the heap to the jitted code, which is then responsible for
305// the deletion of the weak_ptr object.
306template <typename T>
307std::weak_ptr<T> *MakeWeakOnHeap(const std::shared_ptr<T> &shPtr)
308{
309 return new std::weak_ptr<T>(shPtr);
310}
311
312// Same as MakeWeakOnHeap, but create a shared_ptr that makes sure the object is definitely kept alive.
313template <typename T>
314std::shared_ptr<T> *MakeSharedOnHeap(const std::shared_ptr<T> &shPtr)
315{
316 return new std::shared_ptr<T>(shPtr);
317}
318
319bool AtLeastOneEmptyString(const std::vector<std::string_view> strings);
320
321/// Take a shared_ptr<AnyNodeType> and return a shared_ptr<RNodeBase>.
322/// This works for RLoopManager nodes as well as filters and ranges.
323std::shared_ptr<RNodeBase> UpcastNode(std::shared_ptr<RNodeBase> ptr);
324
325ColumnNames_t GetValidatedColumnNames(RLoopManager &lm, const unsigned int nColumns, const ColumnNames_t &columns,
326 const ColumnNames_t &validDefines, RDataSource *ds);
327
328std::vector<std::string> GetValidatedArgTypes(const ColumnNames_t &colNames, const RBookedDefines &defines,
329 TTree *tree, RDataSource *ds, const std::string &context,
330 bool vector2rvec);
331
332std::vector<bool> FindUndefinedDSColumns(const ColumnNames_t &requestedCols, const ColumnNames_t &definedDSCols);
333
335
336template <typename T>
337void AddDSColumnsHelper(const std::string &colName, RLoopManager &lm, RDataSource &ds)
338{
339 if (!ds.HasColumn(colName) || lm.HasDSValuePtrs(colName))
340 return;
341
342 const auto valuePtrs = ds.GetColumnReaders<T>(colName);
343 if (!valuePtrs.empty()) {
344 // we are using the old GetColumnReaders mechanism
345 std::vector<void*> typeErasedValuePtrs(valuePtrs.begin(), valuePtrs.end());
346 lm.AddDSValuePtrs(colName, std::move(typeErasedValuePtrs));
347 }
348}
349
350/// Take list of column names that must be defined, current map of custom columns, current list of defined column names,
351/// and return a new map of custom columns (with the new datasource columns added to it)
352template <typename... ColumnTypes>
353void AddDSColumns(const std::vector<std::string> &requiredCols, RLoopManager &lm, RDataSource &ds,
355{
356 // hack to expand a template parameter pack without c++17 fold expressions.
357 using expander = int[];
358 int i = 0;
359 (void)expander{(AddDSColumnsHelper<ColumnTypes>(requiredCols[i], lm, ds), ++i)..., 0};
360}
361
362// this function is meant to be called by the jitted code generated by BookFilterJit
363template <typename F, typename PrevNode>
364void JitFilterHelper(F &&f, const char **colsPtr, std::size_t colsSize, std::string_view name,
365 std::weak_ptr<RJittedFilter> *wkJittedFilter, std::shared_ptr<PrevNode> *prevNodeOnHeap,
366 RDFInternal::RBookedDefines *defines) noexcept
367{
368 if (wkJittedFilter->expired()) {
369 // The branch of the computation graph that needed this jitted code went out of scope between the type
370 // jitting was booked and the time jitting actually happened. Nothing to do other than cleaning up.
371 delete wkJittedFilter;
372 // defines must be deleted before prevNodeOnHeap because their dtor needs the RLoopManager to be alive
373 // and prevNodeOnHeap is what keeps it alive if the rest of the computation graph is already out of scope
374 delete defines;
375 delete prevNodeOnHeap;
376 return;
377 }
378
379 const ColumnNames_t cols(colsPtr, colsPtr + colsSize);
380 delete[] colsPtr;
381
382 const auto jittedFilter = wkJittedFilter->lock();
383
384 // mock Filter logic -- validity checks and Define-ition of RDataSource columns
385 using Callable_t = typename std::decay<F>::type;
387 using ColTypes_t = typename TTraits::CallableTraits<Callable_t>::arg_types;
388 constexpr auto nColumns = ColTypes_t::list_size;
389 RDFInternal::CheckFilter(f);
390
391 auto &lm = *jittedFilter->GetLoopManagerUnchecked(); // RLoopManager must exist at this time
392 auto ds = lm.GetDataSource();
393
394 if (ds != nullptr)
395 RDFInternal::AddDSColumns(cols, lm, *ds, ColTypes_t());
396
397 jittedFilter->SetFilter(
398 std::unique_ptr<RFilterBase>(new F_t(std::forward<F>(f), cols, *prevNodeOnHeap, *defines, name)));
399 // defines points to the columns structure in the heap, created before the jitted call so that the jitter can
400 // share data after it has lazily compiled the code. Here the data has been used and the memory can be freed.
401 delete defines;
402 delete prevNodeOnHeap;
403 delete wkJittedFilter;
404}
405
406template <typename F>
407void JitDefineHelper(F &&f, const char **colsPtr, std::size_t colsSize, std::string_view name, RLoopManager *lm,
408 std::weak_ptr<RJittedDefine> *wkJittedDefine,
409 RDFInternal::RBookedDefines *defines, std::shared_ptr<RNodeBase> *prevNodeOnHeap) noexcept
410{
411 if (wkJittedDefine->expired()) {
412 // The branch of the computation graph that needed this jitted code went out of scope between the type
413 // jitting was booked and the time jitting actually happened. Nothing to do other than cleaning up.
414 delete wkJittedDefine;
415 // defines must be deleted before prevNodeOnHeap because their dtor needs the RLoopManager to be alive
416 // and prevNodeOnHeap is what keeps it alive if the rest of the computation graph is already out of scope
417 delete defines;
418 delete prevNodeOnHeap;
419 return;
420 }
421
422 const ColumnNames_t cols(colsPtr, colsPtr + colsSize);
423 delete[] colsPtr;
424
425 auto jittedDefine = wkJittedDefine->lock();
426
427 using Callable_t = typename std::decay<F>::type;
429 using ColTypes_t = typename TTraits::CallableTraits<Callable_t>::arg_types;
430 constexpr auto nColumns = ColTypes_t::list_size;
431
432 auto ds = lm->GetDataSource();
433 if (ds != nullptr)
434 RDFInternal::AddDSColumns(cols, *lm, *ds, ColTypes_t());
435
436 // will never actually be used (trumped by jittedDefine->GetTypeName()), but we set it to something meaningful
437 // to help devs debugging
438 const auto dummyType = "jittedCol_t";
439 // use unique_ptr<RDefineBase> instead of make_unique<NewCol_t> to reduce jit/compile-times
440 jittedDefine->SetDefine(std::unique_ptr<RDefineBase>(
441 new NewCol_t(name, dummyType, std::forward<F>(f), cols, lm->GetNSlots(), *defines, lm->GetDSValuePtrs(), ds)));
442
443 // defines points to the columns structure in the heap, created before the jitted call so that the jitter can
444 // share data after it has lazily compiled the code. Here the data has been used and the memory can be freed.
445 delete defines;
446 // prevNodeOnHeap only serves the purpose of keeping the RLoopManager alive so it can be accessed by
447 // defines' destructor in case the rest of the computation graph is gone. Can be safely deleted here.
448 delete prevNodeOnHeap;
449 delete wkJittedDefine;
450}
451
452/// Convenience function invoked by jitted code to build action nodes at runtime
453template <typename ActionTag, typename... ColTypes, typename PrevNodeType, typename HelperArgType>
454void CallBuildAction(std::shared_ptr<PrevNodeType> *prevNodeOnHeap, const char **colsPtr, std::size_t colsSize,
455 const unsigned int nSlots, std::shared_ptr<HelperArgType> *helperArgOnHeap,
456 std::weak_ptr<RJittedAction> *wkJittedActionOnHeap, RDFInternal::RBookedDefines *defines) noexcept
457{
458 if (wkJittedActionOnHeap->expired()) {
459 delete helperArgOnHeap;
460 delete wkJittedActionOnHeap;
461 // defines must be deleted before prevNodeOnHeap because their dtor needs the RLoopManager to be alive
462 // and prevNodeOnHeap is what keeps it alive if the rest of the computation graph is already out of scope
463 delete defines;
464 delete prevNodeOnHeap;
465 return;
466 }
467
468 const ColumnNames_t cols(colsPtr, colsPtr + colsSize);
469 delete[] colsPtr;
470
471 auto jittedActionOnHeap = wkJittedActionOnHeap->lock();
472
473 // if we are here it means we are jitting, if we are jitting the loop manager must be alive
474 auto &prevNodePtr = *prevNodeOnHeap;
475 auto &loopManager = *prevNodePtr->GetLoopManagerUnchecked();
476 using ColTypes_t = TypeList<ColTypes...>;
477 constexpr auto nColumns = ColTypes_t::list_size;
478 auto ds = loopManager.GetDataSource();
479 if (ds != nullptr)
480 RDFInternal::AddDSColumns(cols, loopManager, *ds, ColTypes_t());
481
482 auto actionPtr =
483 BuildAction<ColTypes...>(cols, std::move(*helperArgOnHeap), nSlots, std::move(prevNodePtr), ActionTag{}, *defines);
484 loopManager.AddDataBlockCallback(actionPtr->GetDataBlockCallback());
485 jittedActionOnHeap->SetAction(std::move(actionPtr));
486
487 // defines points to the columns structure in the heap, created before the jitted call so that the jitter can
488 // share data after it has lazily compiled the code. Here the data has been used and the memory can be freed.
489 delete defines;
490
491 delete helperArgOnHeap;
492 delete prevNodeOnHeap;
493 delete wkJittedActionOnHeap;
494}
495
496/// The contained `type` alias is `double` if `T == RInferredType`, `U` if `T == std::container<U>`, `T` otherwise.
497template <typename T, bool Container = RDFInternal::IsDataContainer<T>::value && !std::is_same<T, std::string>::value>
498struct RMinReturnType {
499 using type = T;
500};
501
502template <>
503struct RMinReturnType<RInferredType, false> {
504 using type = double;
505};
506
507template <typename T>
508struct RMinReturnType<T, true> {
509 using type = TTraits::TakeFirstParameter_t<T>;
510};
511
512// return wrapper around f that prepends an `unsigned int slot` parameter
513template <typename R, typename F, typename... Args>
514std::function<R(unsigned int, Args...)> AddSlotParameter(F &f, TypeList<Args...>)
515{
516 return [f](unsigned int, Args... a) mutable -> R { return f(a...); };
517}
518
519template <typename ColType, typename... Rest>
520struct RNeedJittingHelper {
521 static constexpr bool value = RNeedJittingHelper<Rest...>::value;
522};
523
524template <typename... Rest>
525struct RNeedJittingHelper<RInferredType, Rest...> {
526 static constexpr bool value = true;
527};
528
529template <typename T>
530struct RNeedJittingHelper<T> {
531 static constexpr bool value = false;
532};
533
534template <>
535struct RNeedJittingHelper<RInferredType> {
536 static constexpr bool value = true;
537};
538
539template <typename ...ColTypes>
540struct RNeedJitting {
541 static constexpr bool value = RNeedJittingHelper<ColTypes...>::value;
542};
543
544template <>
545struct RNeedJitting<> {
546 static constexpr bool value = false;
547};
548
549///////////////////////////////////////////////////////////////////////////////
550/// Check preconditions for RInterface::Aggregate:
551/// - the aggregator callable must have signature `U(U,T)` or `void(U&,T)`.
552/// - the merge callable must have signature `U(U,U)` or `void(std::vector<U>&)`
553template <typename R, typename Merge, typename U, typename T, typename decayedU = typename std::decay<U>::type,
554 typename mergeArgsNoDecay_t = typename CallableTraits<Merge>::arg_types_nodecay,
555 typename mergeArgs_t = typename CallableTraits<Merge>::arg_types,
556 typename mergeRet_t = typename CallableTraits<Merge>::ret_type>
557void CheckAggregate(TypeList<U, T>)
558{
559 constexpr bool isAggregatorOk =
560 (std::is_same<R, decayedU>::value) || (std::is_same<R, void>::value && std::is_lvalue_reference<U>::value);
561 static_assert(isAggregatorOk, "aggregator function must have signature `U(U,T)` or `void(U&,T)`");
562 constexpr bool isMergeOk =
563 (std::is_same<TypeList<decayedU, decayedU>, mergeArgs_t>::value && std::is_same<decayedU, mergeRet_t>::value) ||
564 (std::is_same<TypeList<std::vector<decayedU> &>, mergeArgsNoDecay_t>::value &&
565 std::is_same<void, mergeRet_t>::value);
566 static_assert(isMergeOk, "merge function must have signature `U(U,U)` or `void(std::vector<U>&)`");
567}
568
569///////////////////////////////////////////////////////////////////////////////
570/// This overload of CheckAggregate is called when the aggregator takes more than two arguments
571template <typename R, typename T>
572void CheckAggregate(T)
573{
574 static_assert(sizeof(T) == 0, "aggregator function must take exactly two arguments");
575}
576
577///////////////////////////////////////////////////////////////////////////////
578/// Check as many template parameters were passed as the number of column names, throw if this is not the case.
579void CheckTypesAndPars(unsigned int nTemplateParams, unsigned int nColumnNames);
580
581/// Return local BranchNames or default BranchNames according to which one should be used
582const ColumnNames_t SelectColumns(unsigned int nArgs, const ColumnNames_t &bl, const ColumnNames_t &defBl);
583
584/// Check whether column names refer to a valid branch of a TTree or have been `Define`d. Return invalid column names.
585ColumnNames_t FindUnknownColumns(const ColumnNames_t &requiredCols, const ColumnNames_t &datasetColumns,
586 const ColumnNames_t &definedCols, const ColumnNames_t &dataSourceColumns);
587
588/// Returns the list of Filters defined in the whole graph
589std::vector<std::string> GetFilterNames(const std::shared_ptr<RLoopManager> &loopManager);
590
591/// Returns the list of Filters defined in the branch
592template <typename NodeType>
593std::vector<std::string> GetFilterNames(const std::shared_ptr<NodeType> &node)
594{
595 std::vector<std::string> filterNames;
596 node->AddFilterName(filterNames);
597 return filterNames;
598}
599
600struct ParsedTreePath {
601 std::string fTreeName;
602 std::string fDirName;
603};
604
605ParsedTreePath ParseTreePath(std::string_view fullTreeName);
606
607// Check if a condition is true for all types
608template <bool...>
609struct TBoolPack;
610
611template <bool... bs>
612using IsTrueForAllImpl_t = typename std::is_same<TBoolPack<bs..., true>, TBoolPack<true, bs...>>;
613
614template <bool... Conditions>
615struct TEvalAnd {
616 static constexpr bool value = IsTrueForAllImpl_t<Conditions...>::value;
617};
618
619// Check if a class is a specialisation of stl containers templates
620// clang-format off
621
622template <typename>
623struct IsList_t : std::false_type {};
624
625template <typename T>
626struct IsList_t<std::list<T>> : std::true_type {};
627
628template <typename>
629struct IsDeque_t : std::false_type {};
630
631template <typename T>
632struct IsDeque_t<std::deque<T>> : std::true_type {};
633// clang-format on
634
635} // namespace RDF
636} // namespace Internal
637
638namespace Detail {
639namespace RDF {
640
641/// The aliased type is `double` if `T == RInferredType`, `U` if `T == container<U>`, `T` otherwise.
642template <typename T>
643using MinReturnType_t = typename RDFInternal::RMinReturnType<T>::type;
644
645template <typename T>
646using MaxReturnType_t = MinReturnType_t<T>;
647
648template <typename T>
649using SumReturnType_t = MinReturnType_t<T>;
650
651} // namespace RDF
652} // namespace Detail
653} // namespace ROOT
654
655/// \endcond
656
657#endif
double
#define d(i)
Definition RSha256.hxx:102
#define f(i)
Definition RSha256.hxx:104
#define g(i)
Definition RSha256.hxx:105
#define a(i)
Definition RSha256.hxx:99
#define h(i)
Definition RSha256.hxx:106
R__EXTERN Int_t gErrorIgnoreLevel
Definition TError.h:129
char name[80]
Definition TGX11.cxx:110
int type
Definition TGX11.cxx:121
typedef void((*Func_t)())
The head node of a RDF computation graph.
void AddDSValuePtrs(const std::string &col, const std::vector< void * > ptrs)
const std::map< std::string, std::vector< void * > > & GetDSValuePtrs() const
RDataSource * GetDataSource() const
bool HasDSValuePtrs(const std::string &col) const
RLoopManager * GetLoopManagerUnchecked() final
Encapsulates the columns defined by the user.
RDataSource defines an API that RDataFrame can use to read arbitrary data formats.
The public interface to the RDataFrame federation of classes.
@ kAllAxes
Definition TH1.h:75
An array of TObjects.
Definition TObjArray.h:37
A TTree represents a columnar dataset.
Definition TTree.h:79
#define F(x, y, z)
std::vector< std::string > ColumnNames_t
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.
ParsedTreePath ParseTreePath(std::string_view fullTreeName)
std::shared_ptr< RNodeBase > UpcastNode(std::shared_ptr< RNodeBase > ptr)
std::vector< std::string > GetFilterNames(const std::shared_ptr< RLoopManager > &loopManager)
std::string PrettyPrintAddr(const void *const addr)
ColumnNames_t GetTopLevelBranchNames(TTree &t)
Get all the top-level branches names, including the ones of the friend trees.
void CheckTypesAndPars(unsigned int nTemplateParams, unsigned int nColumnNames)
std::string DemangleTypeIdName(const std::type_info &typeInfo)
bool AtLeastOneEmptyString(const std::vector< std::string_view > strings)
std::string JitBuildAction(const ColumnNames_t &cols, std::shared_ptr< RDFDetail::RNodeBase > *prevNode, const std::type_info &helperArgType, const std::type_info &at, void *helperArgOnHeap, TTree *tree, const unsigned int nSlots, const RDFInternal::RBookedDefines &customCols, RDataSource *ds, std::weak_ptr< RJittedAction > *jittedActionOnHeap)
ColumnNames_t GetValidatedColumnNames(RLoopManager &lm, const unsigned int nColumns, const ColumnNames_t &columns, const ColumnNames_t &validDefines, RDataSource *ds)
Given the desired number of columns and the user-provided list of columns:
ColumnNames_t FindUnknownColumns(const ColumnNames_t &requiredCols, const ColumnNames_t &datasetColumns, const ColumnNames_t &definedCols, const ColumnNames_t &dataSourceColumns)
void CheckDefine(std::string_view definedCol, TTree *treePtr, const ColumnNames_t &customCols, const std::map< std::string, std::string > &aliasMap, const ColumnNames_t &dataSourceColumns)
ColumnNames_t ConvertRegexToColumns(const ColumnNames_t &colNames, std::string_view columnNameRegexp, std::string_view callerName)
std::vector< std::string > GetValidatedArgTypes(const ColumnNames_t &colNames, const RBookedDefines &defines, TTree *tree, RDataSource *ds, const std::string &context, bool vector2rvec)
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 ...
void BookFilterJit(const std::shared_ptr< RJittedFilter > &jittedFilter, std::shared_ptr< RDFDetail::RNodeBase > *prevNodeOnHeap, std::string_view name, std::string_view expression, const std::map< std::string, std::string > &aliasMap, const ColumnNames_t &branches, const RDFInternal::RBookedDefines &customCols, TTree *tree, RDataSource *ds)
std::shared_ptr< RJittedDefine > BookDefineJit(std::string_view name, std::string_view expression, RLoopManager &lm, RDataSource *ds, const RDFInternal::RBookedDefines &customCols, const ColumnNames_t &branches, std::shared_ptr< RNodeBase > *upcastNodeOnHeap)
double T(double x)
ROOT type_traits extensions.
T Sum(const RVec< T > &v)
Sum elements of an RVec.
Definition RVec.hxx:791
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...
Bool_t IsImplicitMTEnabled()
Returns true if the implicit multi-threading in ROOT is enabled.
Definition TROOT.cxx:556
Short_t Max(Short_t a, Short_t b)
Definition TMathBase.h:212
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:1073
Short_t Min(Short_t a, Short_t b)
Definition TMathBase.h:180
Double_t StdDev(Long64_t n, const T *a, const Double_t *w=0)
Definition TMath.h:530
Definition tree.py:1
A collection of options to steer the creation of the dataset on file.
Lightweight storage for a collection of types.