11#ifndef ROOT_RDF_TINTERFACE
12#define ROOT_RDF_TINTERFACE
38#include <initializer_list>
72template <
typename Proxied,
typename DataSource>
75using RNode = RInterface<::ROOT::Detail::RDF::RNodeBase, void>;
89template <
typename Proxied,
typename DataSource =
void>
99 template <
typename T,
typename W>
126 template <typename T = Proxied, typename std::enable_if<std::is_same<T, RLoopManager>::value,
int>
::type = 0>
186 template <typename F, typename std::enable_if<!std::is_convertible<F, std::string>::value,
int>
::type = 0>
190 RDFInternal::CheckFilter(
f);
191 using ColTypes_t =
typename TTraits::CallableTraits<F>::arg_types;
192 constexpr auto nColumns = ColTypes_t::list_size;
211 template <typename F, typename std::enable_if<!std::is_convertible<F, std::string>::value,
int>
::type = 0>
227 template <
typename F>
252 using BaseNodeType_t =
typename std::remove_pointer<
decltype(upcastNodeOnHeap)>::type::element_type;
254 const auto jittedFilter = std::make_shared<RDFDetail::RJittedFilter>(
fLoopManager,
name);
293 template <typename F, typename std::enable_if<!std::is_convertible<F, std::string>::value,
int>
::type = 0>
296 return DefineImpl<F, RDFDetail::CustomColExtraArgs::None>(
name, std::move(expression), columns);
322 template <
typename F>
325 return DefineImpl<F, RDFDetail::CustomColExtraArgs::Slot>(
name, std::move(expression), columns);
352 template <
typename F>
355 return DefineImpl<F, RDFDetail::CustomColExtraArgs::SlotAndEntry>(
name, std::move(expression), columns);
463 template <
typename... ColumnTypes>
468 return SnapshotImpl<ColumnTypes...>(treename, filename, columnList, options);
489 const auto fullTreeName = treename;
491 treename = parsedTreePath.fTreeName;
492 const auto &dirname = parsedTreePath.fDirName;
494 auto snapHelperArgs = std::make_shared<RDFInternal::SnapshotHelperArgs>(RDFInternal::SnapshotHelperArgs{
495 std::string(filename), std::string(dirname), std::string(treename), columnList, options});
498 auto newRDF = std::make_shared<ROOT::RDataFrame>(fullTreeName, filename, validCols);
500 auto resPtr = CreateAction<RDFInternal::ActionTags::Snapshot, RDFDetail::RInferredType>(
501 validCols, newRDF, snapHelperArgs, validCols.size());
522 std::string_view columnNameRegexp =
"",
530 columnNames.reserve(definedColumns.size() + treeBranchNames.size() + dsColumns.size());
531 columnNames.insert(columnNames.end(), definedColumns.begin(), definedColumns.end());
532 columnNames.insert(columnNames.end(), treeBranchNames.begin(), treeBranchNames.end());
533 columnNames.insert(columnNames.end(), dsColumns.begin(), dsColumns.end());
535 return Snapshot(treename, filename, selectedColumns, options);
553 std::initializer_list<std::string> columnList,
557 return Snapshot(treename, filename, selectedColumns, options);
590 template <
typename... ColumnTypes>
593 auto staticSeq = std::make_index_sequence<
sizeof...(ColumnTypes)>();
594 return CacheImpl<ColumnTypes...>(columnList, staticSeq);
607 if (columnList.empty()) {
608 auto nEntries = *this->
Count();
613 std::stringstream cacheCall;
620 cacheCall <<
"*reinterpret_cast<ROOT::RDF::RInterface<ROOT::Detail::RDF::RLoopManager>*>("
622 <<
") = reinterpret_cast<ROOT::RDF::RInterface<ROOT::Detail::RDF::RNodeBase>*>("
628 for (
const auto &colType : colTypes)
629 cacheCall << colType <<
", ";
630 if (!columnList.empty())
631 cacheCall.seekp(-2, cacheCall.cur);
632 cacheCall <<
">(*reinterpret_cast<std::vector<std::string>*>("
656 columnNames.reserve(definedColumns.size() + treeBranchNames.size() + dsColumns.size());
657 columnNames.insert(columnNames.end(), definedColumns.begin(), definedColumns.end());
658 columnNames.insert(columnNames.end(), treeBranchNames.begin(), treeBranchNames.end());
659 columnNames.insert(columnNames.end(), dsColumns.begin(), dsColumns.end());
661 return Cache(selectedColumns);
673 return Cache(selectedColumns);
697 if (stride == 0 || (end != 0 && end < begin))
698 throw std::runtime_error(
"Range: stride must be strictly greater than 0 and end must be greater than begin.");
702 auto rangePtr = std::make_shared<Range_t>(begin, end, stride,
fProxiedPtr);
735 template <
typename F>
738 using arg_types =
typename TTraits::CallableTraits<
decltype(
f)>::arg_types_nodecay;
739 using ret_type =
typename TTraits::CallableTraits<
decltype(
f)>::ret_type;
740 ForeachSlot(RDFInternal::AddSlotParameter<ret_type>(
f, arg_types()), columns);
765 template <
typename F>
769 constexpr auto nColumns = ColTypes_t::list_size;
774 using Helper_t = RDFInternal::ForeachSlotHelper<F>;
777 auto action = std::make_unique<Action_t>(Helper_t(std::move(
f)), validColumnNames,
fProxiedPtr,
fDefines);
813 template <typename F, typename T = typename TTraits::CallableTraits<F>::ret_type>
817 std::is_default_constructible<T>::value,
818 "reduce object cannot be default-constructed. Please provide an initialisation value (redIdentity)");
819 return Reduce(std::move(
f), columnName, T());
836 template <typename F, typename T = typename TTraits::CallableTraits<F>::ret_type>
858 auto cSPtr = std::make_shared<ULong64_t>(0);
859 using Helper_t = RDFInternal::CountHelper;
864 return MakeResultPtr(cSPtr, *
fLoopManager, std::move(action));
887 template <
typename T,
typename COLL = std::vector<T>>
895 using Helper_t = RDFInternal::TakeHelper<T, T, COLL>;
897 auto valuesPtr = std::make_shared<COLL>();
901 std::make_unique<Action_t>(Helper_t(valuesPtr, nSlots), validColumnNames,
fProxiedPtr,
fDefines);
903 return MakeResultPtr(valuesPtr, *
fLoopManager, std::move(action));
928 template <
typename V = RDFDetail::RInferredType>
935 std::shared_ptr<::TH1D>
h(
nullptr);
937 ROOT::Internal::RDF::RIgnoreErrorLevelRAII iel(
kError);
938 h = model.GetHistogram();
939 h->SetDirectory(
nullptr);
942 if (
h->GetXaxis()->GetXmax() ==
h->GetXaxis()->GetXmin())
943 RDFInternal::HistoUtils<::TH1D>::SetCanExtendAllAxes(*
h);
944 return CreateAction<RDFInternal::ActionTags::Histo1D, V>(validatedColumns,
h,
h);
965 template <
typename V = RDFDetail::RInferredType>
968 const auto h_name = std::string(vName);
969 const auto h_title = h_name +
";" + h_name +
";count";
970 return Histo1D<V>({h_name.c_str(), h_title.c_str(), 128u, 0., 0.}, vName);
992 template <
typename V = RDFDetail::RInferredType,
typename W = RDFDetail::RInferredType>
995 const std::vector<std::string_view> columnViews = {vName, wName};
999 std::shared_ptr<::TH1D>
h(
nullptr);
1001 ROOT::Internal::RDF::RIgnoreErrorLevelRAII iel(
kError);
1004 return CreateAction<RDFInternal::ActionTags::Histo1D, V, W>(userColumns,
h,
h);
1027 template <
typename V = RDFDetail::RInferredType,
typename W = RDFDetail::RInferredType>
1031 std::string str_vName{vName};
1032 std::string str_wName{wName};
1033 const auto h_name = str_vName +
"_weighted_" + str_wName;
1034 const auto h_title = str_vName +
", weights: " + str_wName +
";" + str_vName +
";count * " + str_wName;
1035 return Histo1D<V, W>({h_name.c_str(), h_title.c_str(), 128u, 0., 0.}, vName, wName);
1047 template <
typename V,
typename W>
1050 return Histo1D<V, W>(model,
"",
"");
1077 template <
typename V1 = RDFDetail::RInferredType,
typename V2 = RDFDetail::RInferredType>
1080 std::shared_ptr<::TH2D>
h(
nullptr);
1082 ROOT::Internal::RDF::RIgnoreErrorLevelRAII iel(
kError);
1085 if (!RDFInternal::HistoUtils<::TH2D>::HasAxisLimits(*
h)) {
1086 throw std::runtime_error(
"2D histograms with no axes limits are not supported yet.");
1088 const std::vector<std::string_view> columnViews = {v1Name, v2Name};
1092 return CreateAction<RDFInternal::ActionTags::Histo2D, V1, V2>(userColumns,
h,
h);
1118 template <
typename V1 = RDFDetail::RInferredType,
typename V2 = RDFDetail::RInferredType,
1119 typename W = RDFDetail::RInferredType>
1121 Histo2D(
const TH2DModel &model, std::string_view v1Name, std::string_view v2Name, std::string_view wName)
1123 std::shared_ptr<::TH2D>
h(
nullptr);
1125 ROOT::Internal::RDF::RIgnoreErrorLevelRAII iel(
kError);
1128 if (!RDFInternal::HistoUtils<::TH2D>::HasAxisLimits(*
h)) {
1129 throw std::runtime_error(
"2D histograms with no axes limits are not supported yet.");
1131 const std::vector<std::string_view> columnViews = {v1Name, v2Name, wName};
1135 return CreateAction<RDFInternal::ActionTags::Histo2D, V1, V2, W>(userColumns,
h,
h);
1138 template <
typename V1,
typename V2,
typename W>
1141 return Histo2D<V1, V2, W>(model,
"",
"",
"");
1168 template <
typename V1 = RDFDetail::RInferredType,
typename V2 = RDFDetail::RInferredType,
1169 typename V3 = RDFDetail::RInferredType>
1171 std::string_view v3Name =
"")
1173 std::shared_ptr<::TH3D>
h(
nullptr);
1175 ROOT::Internal::RDF::RIgnoreErrorLevelRAII iel(
kError);
1178 if (!RDFInternal::HistoUtils<::TH3D>::HasAxisLimits(*
h)) {
1179 throw std::runtime_error(
"3D histograms with no axes limits are not supported yet.");
1181 const std::vector<std::string_view> columnViews = {v1Name, v2Name, v3Name};
1185 return CreateAction<RDFInternal::ActionTags::Histo3D, V1, V2, V3>(userColumns,
h,
h);
1215 template <
typename V1 = RDFDetail::RInferredType,
typename V2 = RDFDetail::RInferredType,
1216 typename V3 = RDFDetail::RInferredType,
typename W = RDFDetail::RInferredType>
1218 std::string_view v3Name, std::string_view wName)
1220 std::shared_ptr<::TH3D>
h(
nullptr);
1222 ROOT::Internal::RDF::RIgnoreErrorLevelRAII iel(
kError);
1225 if (!RDFInternal::HistoUtils<::TH3D>::HasAxisLimits(*
h)) {
1226 throw std::runtime_error(
"3D histograms with no axes limits are not supported yet.");
1228 const std::vector<std::string_view> columnViews = {v1Name, v2Name, v3Name, wName};
1232 return CreateAction<RDFInternal::ActionTags::Histo3D, V1, V2, V3, W>(userColumns,
h,
h);
1235 template <
typename V1,
typename V2,
typename V3,
typename W>
1238 return Histo3D<V1, V2, V3, W>(model,
"",
"",
"",
"");
1266 template <
typename V1 = RDFDetail::RInferredType,
typename V2 = RDFDetail::RInferredType>
1269 auto graph = std::make_shared<::TGraph>();
1270 const std::vector<std::string_view> columnViews = {v1Name, v2Name};
1278 if (!(validatedColumns[0].empty() && validatedColumns[1].empty())) {
1279 const auto g_name = std::string(v1Name) +
"_vs_" + std::string(v2Name);
1280 const auto g_title = std::string(v1Name) +
" vs " + std::string(v2Name);
1281 graph->SetNameTitle(g_name.c_str(), g_title.c_str());
1282 graph->GetXaxis()->SetTitle(std::string(v1Name).c_str());
1283 graph->GetYaxis()->SetTitle(std::string(v2Name).c_str());
1286 return CreateAction<RDFInternal::ActionTags::Graph, V1, V2>(validatedColumns,
graph,
graph);
1309 template <
typename V1 = RDFDetail::RInferredType,
typename V2 = RDFDetail::RInferredType>
1313 std::shared_ptr<::TProfile>
h(
nullptr);
1315 ROOT::Internal::RDF::RIgnoreErrorLevelRAII iel(
kError);
1319 if (!RDFInternal::HistoUtils<::TProfile>::HasAxisLimits(*
h)) {
1320 throw std::runtime_error(
"Profiles with no axes limits are not supported yet.");
1322 const std::vector<std::string_view> columnViews = {v1Name, v2Name};
1326 return CreateAction<RDFInternal::ActionTags::Profile1D, V1, V2>(userColumns,
h,
h);
1352 template <
typename V1 = RDFDetail::RInferredType,
typename V2 = RDFDetail::RInferredType,
1353 typename W = RDFDetail::RInferredType>
1357 std::shared_ptr<::TProfile>
h(
nullptr);
1359 ROOT::Internal::RDF::RIgnoreErrorLevelRAII iel(
kError);
1363 if (!RDFInternal::HistoUtils<::TProfile>::HasAxisLimits(*
h)) {
1364 throw std::runtime_error(
"Profile histograms with no axes limits are not supported yet.");
1366 const std::vector<std::string_view> columnViews = {v1Name, v2Name, wName};
1370 return CreateAction<RDFInternal::ActionTags::Profile1D, V1, V2, W>(userColumns,
h,
h);
1373 template <
typename V1,
typename V2,
typename W>
1376 return Profile1D<V1, V2, W>(model,
"",
"",
"");
1403 template <
typename V1 = RDFDetail::RInferredType,
typename V2 = RDFDetail::RInferredType,
1404 typename V3 = RDFDetail::RInferredType>
1406 std::string_view v2Name =
"", std::string_view v3Name =
"")
1408 std::shared_ptr<::TProfile2D>
h(
nullptr);
1410 ROOT::Internal::RDF::RIgnoreErrorLevelRAII iel(
kError);
1414 if (!RDFInternal::HistoUtils<::TProfile2D>::HasAxisLimits(*
h)) {
1415 throw std::runtime_error(
"2D profiles with no axes limits are not supported yet.");
1417 const std::vector<std::string_view> columnViews = {v1Name, v2Name, v3Name};
1421 return CreateAction<RDFInternal::ActionTags::Profile2D, V1, V2, V3>(userColumns,
h,
h);
1450 template <
typename V1 = RDFDetail::RInferredType,
typename V2 = RDFDetail::RInferredType,
1451 typename V3 = RDFDetail::RInferredType,
typename W = RDFDetail::RInferredType>
1453 std::string_view v3Name, std::string_view wName)
1455 std::shared_ptr<::TProfile2D>
h(
nullptr);
1457 ROOT::Internal::RDF::RIgnoreErrorLevelRAII iel(
kError);
1461 if (!RDFInternal::HistoUtils<::TProfile2D>::HasAxisLimits(*
h)) {
1462 throw std::runtime_error(
"2D profiles with no axes limits are not supported yet.");
1464 const std::vector<std::string_view> columnViews = {v1Name, v2Name, v3Name, wName};
1468 return CreateAction<RDFInternal::ActionTags::Profile2D, V1, V2, V3, W>(userColumns,
h,
h);
1471 template <
typename V1,
typename V2,
typename V3,
typename W>
1474 return Profile2D<V1, V2, V3, W>(model,
"",
"",
"",
"");
1501 template <
typename FirstColumn,
typename... OtherColumns,
typename T>
1504 auto h = std::make_shared<T>(std::forward<T>(model));
1505 if (!RDFInternal::HistoUtils<T>::HasAxisLimits(*
h)) {
1506 throw std::runtime_error(
"The absence of axes limits is not supported yet.");
1508 return CreateAction<RDFInternal::ActionTags::Fill, FirstColumn, OtherColumns...>(columnList,
h,
h);
1530 template <
typename T>
1533 auto h = std::make_shared<T>(std::forward<T>(model));
1534 if (!RDFInternal::HistoUtils<T>::HasAxisLimits(*
h)) {
1535 throw std::runtime_error(
"The absence of axes limits is not supported yet.");
1537 return CreateAction<RDFInternal::ActionTags::Fill, RDFDetail::RInferredType>(columnList,
h,
h, columnList.size());
1555 template <
typename V = RDFDetail::RInferredType>
1559 if (!value.empty()) {
1560 columns.emplace_back(std::string(value));
1563 if (std::is_same<V, RDFDetail::RInferredType>::value) {
1566 return Fill<V>(
TStatistic(), validColumnNames);
1587 template <
typename V = RDFDetail::RInferredType,
typename W = RDFDetail::RInferredType>
1590 ColumnNames_t columns{std::string(value), std::string(weight)};
1591 constexpr auto vIsInferred = std::is_same<V, RDFDetail::RInferredType>::value;
1592 constexpr auto wIsInferred = std::is_same<W, RDFDetail::RInferredType>::value;
1598 if (vIsInferred && wIsInferred) {
1600 }
else if (vIsInferred != wIsInferred) {
1601 std::string error(
"The ");
1602 error += vIsInferred ?
"value " :
"weight ";
1603 error +=
"column type is explicit, while the ";
1604 error += vIsInferred ?
"weight " :
"value ";
1605 error +=
" is specified to be inferred. This case is not supported: please specify both types or none.";
1606 throw std::runtime_error(error);
1608 return Fill<V, W>(
TStatistic(), validColumnNames);
1633 template <
typename T = RDFDetail::RInferredType>
1637 using RetType_t = RDFDetail::MinReturnType_t<T>;
1638 auto minV = std::make_shared<RetType_t>(std::numeric_limits<RetType_t>::max());
1639 return CreateAction<RDFInternal::ActionTags::Min, T>(userColumns, minV, minV);
1663 template <
typename T = RDFDetail::RInferredType>
1667 using RetType_t = RDFDetail::MaxReturnType_t<T>;
1668 auto maxV = std::make_shared<RetType_t>(std::numeric_limits<RetType_t>::lowest());
1669 return CreateAction<RDFInternal::ActionTags::Max, T>(userColumns, maxV, maxV);
1692 template <
typename T = RDFDetail::RInferredType>
1696 auto meanV = std::make_shared<double>(0);
1697 return CreateAction<RDFInternal::ActionTags::Mean, T>(userColumns, meanV, meanV);
1720 template <
typename T = RDFDetail::RInferredType>
1724 auto stdDeviationV = std::make_shared<double>(0);
1725 return CreateAction<RDFInternal::ActionTags::StdDev, T>(userColumns, stdDeviationV, stdDeviationV);
1751 template <
typename T = RDFDetail::RInferredType>
1753 Sum(std::string_view columnName =
"",
1754 const RDFDetail::SumReturnType_t<T> &initValue = RDFDetail::SumReturnType_t<T>{})
1757 auto sumV = std::make_shared<RDFDetail::SumReturnType_t<T>>(initValue);
1758 return CreateAction<RDFInternal::ActionTags::Sum, T>(userColumns, sumV, sumV);
1788 bool returnEmptyReport =
false;
1794 if (std::is_same<Proxied, RLoopManager>::value &&
fDefines.
GetNames().size() > 4)
1795 returnEmptyReport =
true;
1797 auto rep = std::make_shared<RCutFlowReport>();
1798 using Helper_t = RDFInternal::ReportHelper<Proxied>;
1805 return MakeResultPtr(rep, *
fLoopManager, std::move(action));
1825 auto addIfNotInternal = [&allColumns](std::string_view colName) {
1827 allColumns.emplace_back(colName);
1832 std::for_each(columnNames.begin(), columnNames.end(), addIfNotInternal);
1837 allColumns.insert(allColumns.end(), branchNames.begin(), branchNames.end());
1842 allColumns.insert(allColumns.end(), dsColNames.begin(), dsColNames.end());
1863 auto col = std::string(column);
1867 const auto it = aliasMap.find(col);
1868 if (it != aliasMap.end())
1873 const bool convertVector2RVec =
true;
1875 convertVector2RVec);
1915 for (
auto column : columns) {
1917 definedColumns.emplace_back(column.first);
1920 return definedColumns;
1944 const auto branchNamesEnd = branchNames.end();
1945 if (branchNamesEnd != std::find(branchNames.begin(), branchNamesEnd, columnName))
2030 template <typename AccFun, typename MergeFun, typename R = typename TTraits::CallableTraits<AccFun>::ret_type,
2031 typename ArgTypes =
typename TTraits::CallableTraits<AccFun>::arg_types,
2032 typename ArgTypesNoDecay =
typename TTraits::CallableTraits<AccFun>::arg_types_nodecay,
2033 typename U = TTraits::TakeFirstParameter_t<ArgTypes>,
2034 typename T = TTraits::TakeFirstParameter_t<TTraits::RemoveFirstParameter_t<ArgTypes>>>
2037 RDFInternal::CheckAggregate<R, MergeFun>(ArgTypesNoDecay());
2043 auto accObjPtr = std::make_shared<U>(aggIdentity);
2044 using Helper_t = RDFInternal::AggregateHelper<AccFun, MergeFun, R, T, U>;
2046 auto action = std::make_unique<Action_t>(
2047 Helper_t(std::move(aggregator), std::move(merger), accObjPtr,
fLoopManager->
GetNSlots()), validColumnNames,
2050 return MakeResultPtr(accObjPtr, *
fLoopManager, std::move(action));
2066 template <typename AccFun, typename MergeFun, typename R = typename TTraits::CallableTraits<AccFun>::ret_type,
2067 typename ArgTypes =
typename TTraits::CallableTraits<AccFun>::arg_types,
2068 typename U = TTraits::TakeFirstParameter_t<ArgTypes>,
2069 typename T = TTraits::TakeFirstParameter_t<TTraits::RemoveFirstParameter_t<ArgTypes>>>
2073 std::is_default_constructible<U>::value,
2074 "aggregated object cannot be default-constructed. Please provide an initialisation value (aggIdentity)");
2075 return Aggregate(std::move(aggregator), std::move(merger), columnName, U());
2112 template <
typename... ColumnTypes,
typename Helper>
2115 constexpr auto nColumns =
sizeof...(ColumnTypes);
2121 using AH = RDFDetail::RActionImpl<Helper>;
2122 static_assert(std::is_base_of<AH, Helper>::value && std::is_convertible<Helper *, AH *>::value,
2123 "Action helper of type T must publicly inherit from ROOT::Detail::RDF::RActionImpl<T>");
2126 auto resPtr = helper.GetResultPtr();
2130 auto action = std::make_unique<Action_t>(Helper(std::forward<Helper>(helper)), validColumnNames,
fProxiedPtr,
2134 return MakeResultPtr(resPtr, *
fLoopManager, std::move(action));
2160 template <
typename... ColumnTypes>
2165 auto displayer = std::make_shared<RDFInternal::RDisplay>(columnList,
GetColumnTypeNamesList(columnList), nRows);
2166 return CreateAction<RDFInternal::ActionTags::Display, ColumnTypes...>(columnList, displayer, displayer);
2180 auto displayer = std::make_shared<RDFInternal::RDisplay>(columnList,
GetColumnTypeNamesList(columnList), nRows);
2181 return CreateAction<RDFInternal::ActionTags::Display, RDFDetail::RInferredType>(columnList, displayer, displayer,
2198 return Display(selectedColumns, nRows);
2211 return Display(selectedColumns, nRows);
2220 const std::string entryColName =
"rdfentry_";
2221 const std::string entryColType =
"ULong64_t";
2222 auto entryColGen = [](
unsigned int,
ULong64_t entry) {
return entry; };
2223 using NewColEntry_t =
2224 RDFDetail::RDefine<
decltype(entryColGen), RDFDetail::CustomColExtraArgs::SlotAndEntry>;
2226 auto entryColumn = std::make_shared<NewColEntry_t>(entryColName, entryColType, std::move(entryColGen),
2229 newCols.
AddColumn(entryColumn, entryColName);
2232 const std::string slotColName =
"rdfslot_";
2233 const std::string slotColType =
"unsigned int";
2234 auto slotColGen = [](
unsigned int slot) {
return slot; };
2235 using NewColSlot_t =
RDFDetail::RDefine<
decltype(slotColGen), RDFDetail::CustomColExtraArgs::Slot>;
2237 auto slotColumn = std::make_shared<NewColSlot_t>(slotColName, slotColType, std::move(slotColGen),
ColumnNames_t{},
2240 newCols.
AddColumn(slotColumn, slotColName);
2252 std::vector<std::string> types;
2254 for (
auto column : columnList) {
2263 std::string error(callerName);
2264 error +=
" was called with ImplicitMT enabled, but multi-thread is not supported.";
2265 throw std::runtime_error(error);
2274 template <
typename ActionTag,
typename... ColTypes,
typename ActionResultType,
2275 typename HelperArgType = ActionResultType,
2276 typename std::enable_if<!RDFInternal::RNeedJitting<ColTypes...>::value,
int>
::type = 0>
2278 const std::shared_ptr<HelperArgType> &helperArg)
2280 constexpr auto nColumns =
sizeof...(ColTypes);
2288 RDFInternal::BuildAction<ColTypes...>(validColumnNames, helperArg, nSlots,
fProxiedPtr, ActionTag{},
fDefines);
2298 template <
typename ActionTag,
typename... ColTypes,
typename ActionResultType,
2299 typename HelperArgType = ActionResultType,
2300 typename std::enable_if<RDFInternal::RNeedJitting<ColTypes...>::value,
int>
::type = 0>
2302 const std::shared_ptr<HelperArgType> &helperArg,
const int nColumns = -1)
2304 auto realNColumns = (nColumns > -1 ? nColumns :
sizeof...(ColTypes));
2310 auto helperArgOnHeap = RDFInternal::MakeSharedOnHeap(helperArg);
2313 using BaseNodeType_t =
typename std::remove_pointer<
decltype(upcastNodeOnHeap)>::type::element_type;
2316 const auto jittedAction = std::make_shared<RDFInternal::RJittedAction>(*
fLoopManager);
2317 auto jittedActionOnHeap = RDFInternal::MakeWeakOnHeap(jittedAction);
2320 validColumnNames, upcastNodeOnHeap,
typeid(std::shared_ptr<HelperArgType>),
typeid(ActionTag), helperArgOnHeap,
2324 return MakeResultPtr(
r, *
fLoopManager, std::move(jittedAction));
2327 template <typename F, typename DefineType, typename RetType = typename TTraits::CallableTraits<F>::ret_type>
2335 using ArgTypes_t =
typename TTraits::CallableTraits<F>::arg_types;
2336 using ColTypesTmp_t =
typename RDFInternal::RemoveFirstParameterIf<
2337 std::is_same<DefineType, RDFDetail::CustomColExtraArgs::Slot>::value, ArgTypes_t>
::type;
2338 using ColTypes_t =
typename RDFInternal::RemoveFirstTwoParametersIf<
2339 std::is_same<DefineType, RDFDetail::CustomColExtraArgs::SlotAndEntry>::value, ColTypesTmp_t>
::type;
2341 constexpr auto nColumns = ColTypes_t::list_size;
2348 if (retTypeName.empty()) {
2352 retTypeName =
"CLING_UNKNOWN_TYPE_" + demangledType;
2357 std::make_shared<NewCol_t>(
name, retTypeName, std::forward<F>(expression), validColumnNames,
2365 return newInterface;
2371 template <typename F, typename DefineType, typename RetType = typename TTraits::CallableTraits<F>::ret_type,
2372 bool IsFStringConv = std::is_convertible<F, std::string>::value,
2373 bool IsRetTypeDefConstr = std::is_default_constructible<RetType>::value>
2374 typename std::enable_if<!IsFStringConv && !IsRetTypeDefConstr, RInterface<Proxied, DS_t>>
::type
2377 static_assert(std::is_default_constructible<typename TTraits::CallableTraits<F>::ret_type>::value,
2378 "Error in `Define`: type returned by expression is not default-constructible");
2382 template <
typename... ColumnTypes>
2392 const auto &treename = parsedTreePath.fTreeName;
2393 const auto &dirname = parsedTreePath.fDirName;
2395 auto snapHelperArgs = std::make_shared<RDFInternal::SnapshotHelperArgs>(RDFInternal::SnapshotHelperArgs{
2396 std::string(filename), std::string(dirname), std::string(treename), columnList, options});
2399 auto newRDF = std::make_shared<ROOT::RDataFrame>(fullTreeName, filename, validCols);
2401 auto resPtr =
CreateAction<RDFInternal::ActionTags::Snapshot, ColumnTypes...>(validCols, newRDF, snapHelperArgs);
2410 template <
typename... ColTypes, std::size_t... S>
2414 constexpr bool areCopyConstructible =
2415 RDFInternal::TEvalAnd<std::is_copy_constructible<ColTypes>::value...>::value;
2416 static_assert(areCopyConstructible,
"Columns of a type which is not copy constructible cannot be cached yet.");
2420 auto colHolders = std::make_tuple(Take<ColTypes>(columnList[S])...);
2421 auto ds = std::make_unique<
RLazyDS<ColTypes...>>(std::make_pair(columnList[S], std::get<S>(colHolders))...);
2447 template <
typename... ColumnTypes>
unsigned long long ULong64_t
The head node of a RDF computation graph.
const std::map< std::string, std::string > & GetAliasMap() const
const ColumnNames_t & GetBranchNames()
Return all valid TTree::Branch names (caching results for subsequent calls).
void ToJitExec(const std::string &) const
unsigned int GetNRuns() const
void Run()
Start the event loop with a different mechanism depending on IMT/no IMT, data source/no data source.
void AddColumnAlias(const std::string &alias, const std::string &colName)
const std::map< std::string, std::vector< void * > > & GetDSValuePtrs() const
RDataSource * GetDataSource() const
unsigned int GetNSlots() const
void Book(RDFInternal::RActionBase *actionPtr)
void Jit()
Add RDF nodes that require just-in-time compilation to the computation graph.
void AddDataBlockCallback(std::function< void(unsigned int)> &&callback)
Helper class that provides the operation graph nodes.
A RDataFrame node that produces a result.
Encapsulates the columns defined by the user.
bool HasName(std::string_view name) const
Check if the provided name is tracked in the names list.
void AddColumn(const std::shared_ptr< RDFDetail::RDefineBase > &column, std::string_view name)
Add a new booked column.
const RDefineBasePtrMap_t & GetColumns() const
Returns the list of the pointers to the defined columns.
ColumnNames_t GetNames() const
Returns the list of the names of the defined columns.
void AddName(std::string_view name)
Add a new name to the list returned by GetNames without booking a new column.
RDataSource defines an API that RDataFrame can use to read arbitrary data formats.
virtual bool HasColumn(std::string_view colName) const =0
Checks if the dataset has a certain column.
virtual const std::vector< std::string > & GetColumnNames() const =0
Returns a reference to the collection of the dataset's column names.
The public interface to the RDataFrame federation of classes.
RInterface(const RInterface &)=default
Copy-ctor for RInterface.
RResultPtr<::TH1D > Histo1D(std::string_view vName, std::string_view wName)
Fill and return a one-dimensional histogram with the weighted values of a column (lazy action)
RInterface(const std::shared_ptr< Proxied > &proxied, RLoopManager &lm, const RDFInternal::RBookedDefines &columns, RDataSource *ds)
RResultPtr<::TH1D > Histo1D(const TH1DModel &model={"", "", 128u, 0., 0.})
Fill and return a one-dimensional histogram with the weighted values of a column (lazy action)
RResultPtr< T > Fill(T &&model, const ColumnNames_t &columnList)
Return an object of type T on which T::Fill will be called once per event (lazy action)
RResultPtr<::TH2D > Histo2D(const TH2DModel &model)
RResultPtr<::TProfile > Profile1D(const TProfile1DModel &model, std::string_view v1Name="", std::string_view v2Name="")
Fill and return a one-dimensional profile (lazy action)
RResultPtr< RInterface< RLoopManager > > Snapshot(std::string_view treename, std::string_view filename, std::string_view columnNameRegexp="", const RSnapshotOptions &options=RSnapshotOptions())
Save selected columns to disk, in a new TTree treename in file filename.
RResultPtr< TStatistic > Stats(std::string_view value="")
Return a TStatistic object, filled once per event (lazy action)
RLoopManager * GetLoopManager() const
RResultPtr<::TGraph > Graph(std::string_view v1Name="", std::string_view v2Name="")
Fill and return a graph (lazy action)
RInterface< Proxied, DS_t > DefineSlot(std::string_view name, F expression, const ColumnNames_t &columns={})
Creates a custom column with a value dependent on the processing slot.
RResultPtr< double > StdDev(std::string_view columnName="")
Return the unbiased standard deviation of processed column values (lazy action)
unsigned int GetNSlots() const
Gets the number of data processing slots.
RInterface(const std::shared_ptr< Proxied > &proxied)
Only enabled when building a RInterface<RLoopManager>
void ForeachSlot(F f, const ColumnNames_t &columns={})
Execute a user-defined function requiring a processing slot index on each entry (instant action)
RResultPtr< ActionResultType > CreateAction(const ColumnNames_t &columns, const std::shared_ptr< ActionResultType > &r, const std::shared_ptr< HelperArgType > &helperArg, const int nColumns=-1)
Create RAction object, return RResultPtr for the action Overload for the case in which one or more co...
RInterface< RLoopManager > Cache(const ColumnNames_t &columnList)
Save selected columns in memory.
RResultPtr< RDisplay > Display(std::initializer_list< std::string > columnList, const int &nRows=5)
Provides a representation of the columns in the dataset.
RInterface< Proxied, DS_t > Define(std::string_view name, F expression, const ColumnNames_t &columns={})
Creates a custom column.
RResultPtr< TStatistic > Stats(std::string_view value, std::string_view weight)
Return a TStatistic object, filled once per event (lazy action)
RDataSource * fDataSource
Non-owning pointer to a data-source object. Null if no data-source. RLoopManager has ownership of the...
RResultPtr<::TH2D > Histo2D(const TH2DModel &model, std::string_view v1Name="", std::string_view v2Name="")
Fill and return a two-dimensional histogram (lazy action)
RResultPtr< RInterface< RLoopManager > > SnapshotImpl(std::string_view fullTreeName, std::string_view filename, const ColumnNames_t &columnList, const RSnapshotOptions &options)
RResultPtr< ActionResultType > CreateAction(const ColumnNames_t &columns, const std::shared_ptr< ActionResultType > &r, const std::shared_ptr< HelperArgType > &helperArg)
Create RAction object, return RResultPtr for the action Overload for the case in which all column typ...
ColumnNames_t GetValidatedColumnNames(const unsigned int nColumns, const ColumnNames_t &columns)
Prepare the call to the GetValidatedColumnNames routine, making sure that GetBranchNames,...
RResultPtr<::TProfile > Profile1D(const TProfile1DModel &model)
RInterface< RDFDetail::RFilter< F, Proxied >, DS_t > Filter(F f, const std::initializer_list< std::string > &columns)
Append a filter to the call graph.
RResultPtr< double > Mean(std::string_view columnName="")
Return the mean of processed column values (lazy action)
RResultPtr< RInterface< RLoopManager > > Snapshot(std::string_view treename, std::string_view filename, std::initializer_list< std::string > columnList, const RSnapshotOptions &options=RSnapshotOptions())
Save selected columns to disk, in a new TTree treename in file filename.
std::enable_if<!IsFStringConv &&!IsRetTypeDefConstr, RInterface< Proxied, DS_t > >::type DefineImpl(std::string_view, F, const ColumnNames_t &)
RInterface< Proxied, DS_t > Alias(std::string_view alias, std::string_view columnName)
Allow to refer to a column with a different name.
RInterface< RLoopManager > Cache(const ColumnNames_t &columnList)
Save selected columns in memory.
RDFInternal::RBookedDefines fDefines
Contains the custom columns defined up to this node.
RInterface< RLoopManager > Cache(std::string_view columnNameRegexp="")
Save selected columns in memory.
RResultPtr< RDisplay > Display(std::string_view columnNameRegexp="", const int &nRows=5)
Provides a representation of the columns in the dataset.
RLoopManager * fLoopManager
friend class RDFInternal::GraphDrawing::GraphCreatorHelper
RResultPtr<::TH2D > Histo2D(const TH2DModel &model, std::string_view v1Name, std::string_view v2Name, std::string_view wName)
Fill and return a weighted two-dimensional histogram (lazy action)
RInterface & operator=(const RInterface &)=default
Copy-assignment operator for RInterface.
RResultPtr< RDFDetail::SumReturnType_t< T > > Sum(std::string_view columnName="", const RDFDetail::SumReturnType_t< T > &initValue=RDFDetail::SumReturnType_t< T >{})
Return the sum of processed column values (lazy action)
RResultPtr< ULong64_t > Count()
Return the number of entries processed (lazy action)
RResultPtr< T > Fill(T &&model, const ColumnNames_t &columnList)
Return an object of type T on which T::Fill will be called once per event (lazy action)
RInterface< Proxied, DS_t > Define(std::string_view name, std::string_view expression)
Creates a custom column.
std::shared_ptr< Proxied > fProxiedPtr
Smart pointer to the graph node encapsulated by this RInterface.
RResultPtr<::TH1D > Histo1D(std::string_view vName)
Fill and return a one-dimensional histogram with the values of a column (lazy action)
ColumnNames_t GetColumnNames()
Returns the names of the available columns.
RResultPtr< RDisplay > Display(const ColumnNames_t &columnList, const int &nRows=5)
Provides a representation of the columns in the dataset.
RResultPtr<::TH1D > Histo1D(const TH1DModel &model, std::string_view vName, std::string_view wName)
Fill and return a one-dimensional histogram with the weighted values of a column (lazy action)
RInterface< RLoopManager > CacheImpl(const ColumnNames_t &columnList, std::index_sequence< S... >)
Implementation of cache.
RInterface< RDFDetail::RRange< Proxied >, DS_t > Range(unsigned int end)
Creates a node that filters entries based on range.
RResultPtr< COLL > Take(std::string_view column="")
Return a collection of values of a column (lazy action, returns a std::vector by default)
RInterface< RLoopManager > Cache(std::initializer_list< std::string > columnList)
Save selected columns in memory.
void CheckAndFillDSColumns(ColumnNames_t validCols, TTraits::TypeList< ColumnTypes... > typeList)
RResultPtr<::TProfile2D > Profile2D(const TProfile2DModel &model, std::string_view v1Name="", std::string_view v2Name="", std::string_view v3Name="")
Fill and return a two-dimensional profile (lazy action)
RResultPtr< typename Helper::Result_t > Book(Helper &&helper, const ColumnNames_t &columns={})
Book execution of a custom action using a user-defined helper object.
RResultPtr< RDisplay > Display(const ColumnNames_t &columnList, const int &nRows=5)
Provides a representation of the columns in the dataset.
const std::shared_ptr< Proxied > & GetProxiedPtr() const
RResultPtr<::TH3D > Histo3D(const TH3DModel &model, std::string_view v1Name="", std::string_view v2Name="", std::string_view v3Name="")
Fill and return a three-dimensional histogram (lazy action)
RResultPtr< RInterface< RLoopManager > > Snapshot(std::string_view treename, std::string_view filename, const ColumnNames_t &columnList, const RSnapshotOptions &options=RSnapshotOptions())
Save selected columns to disk, in a new TTree treename in file filename.
RResultPtr< RCutFlowReport > Report()
Gather filtering statistics.
RResultPtr<::TProfile2D > Profile2D(const TProfile2DModel &model, std::string_view v1Name, std::string_view v2Name, std::string_view v3Name, std::string_view wName)
Fill and return a two-dimensional profile (lazy action)
RResultPtr< RInterface< RLoopManager > > Snapshot(std::string_view treename, std::string_view filename, const ColumnNames_t &columnList, const RSnapshotOptions &options=RSnapshotOptions())
Save selected columns to disk, in a new TTree treename in file filename.
RResultPtr< U > Aggregate(AccFun aggregator, MergeFun merger, std::string_view columnName="")
Execute a user-defined accumulation operation on the processed column values in each processing slot.
RInterface< Proxied, DS_t > DefineSlotEntry(std::string_view name, F expression, const ColumnNames_t &columns={})
Creates a custom column with a value dependent on the processing slot and the current entry.
RResultPtr< RDFDetail::MinReturnType_t< T > > Min(std::string_view columnName="")
Return the minimum of processed column values (lazy action)
RResultPtr< T > Reduce(F f, std::string_view columnName="")
Execute a user-defined reduce operation on the values of a column.
void Foreach(F f, const ColumnNames_t &columns={})
Execute a user-defined function on each entry (instant action)
RInterface< RDFDetail::RJittedFilter, DS_t > Filter(std::string_view expression, std::string_view name="")
Append a filter to the call graph.
RResultPtr<::TProfile2D > Profile2D(const TProfile2DModel &model)
std::enable_if< std::is_default_constructible< RetType >::value, RInterface< Proxied, DS_t > >::type DefineImpl(std::string_view name, F &&expression, const ColumnNames_t &columns)
std::string GetColumnType(std::string_view column)
Return the type of a given column as a string.
ColumnNames_t GetDefinedColumnNames()
Returns the names of the defined columns.
RInterface< RDFDetail::RFilter< F, Proxied >, DS_t > Filter(F f, const ColumnNames_t &columns={}, std::string_view name="")
Append a filter to the call graph.
RResultPtr< U > Aggregate(AccFun aggregator, MergeFun merger, std::string_view columnName, const U &aggIdentity)
Execute a user-defined accumulation operation on the processed column values in each processing slot.
RInterface(RInterface &&)=default
Move-ctor for RInterface.
RResultPtr< T > Reduce(F f, std::string_view columnName, const T &redIdentity)
Execute a user-defined reduce operation on the values of a column.
void CheckIMTDisabled(std::string_view callerName)
unsigned int GetNRuns() const
Gets the number of event loops run.
RResultPtr<::TH3D > Histo3D(const TH3DModel &model, std::string_view v1Name, std::string_view v2Name, std::string_view v3Name, std::string_view wName)
Fill and return a three-dimensional histogram (lazy action)
bool HasColumn(std::string_view columnName)
Checks if a column is present in the dataset.
RDFDetail::ColumnNames_t ColumnNames_t
RInterface< RDFDetail::RFilter< F, Proxied >, DS_t > Filter(F f, std::string_view name)
Append a filter to the call graph.
RInterface< RDFDetail::RRange< Proxied >, DS_t > Range(unsigned int begin, unsigned int end, unsigned int stride=1)
Creates a node that filters entries based on range: [begin, end)
std::vector< std::string > GetColumnTypeNamesList(const ColumnNames_t &columnList)
std::vector< std::string > GetFilterNames()
Returns the names of the filters created.
RResultPtr<::TH1D > Histo1D(const TH1DModel &model={"", "", 128u, 0., 0.}, std::string_view vName="")
Fill and return a one-dimensional histogram with the values of a column (lazy action)
RResultPtr<::TProfile > Profile1D(const TProfile1DModel &model, std::string_view v1Name, std::string_view v2Name, std::string_view wName)
Fill and return a one-dimensional profile (lazy action)
RResultPtr<::TH3D > Histo3D(const TH3DModel &model)
RResultPtr< RDFDetail::MaxReturnType_t< T > > Max(std::string_view columnName="")
Return the maximum of processed column values (lazy action)
A RDataSource implementation which is built on top of result proxies.
Smart pointer for the return type of actions.
ROOT's RDataFrame offers a high level interface for analyses of data stored in TTrees,...
typename RemoveFirstParameter< T >::type RemoveFirstParameter_t
Small helper to keep current directory context.
A TGraph is an object made of two arrays X and Y with npoints each.
Statistical variable, defined by its mean and variance (RMS).
std::vector< std::string > ColumnNames_t
ParsedTreePath ParseTreePath(std::string_view fullTreeName)
std::vector< std::string > GetBranchNames(TTree &t, bool allowDuplicates=true)
Get all the branches names, including the ones of the friend trees.
std::string ColumnName2ColumnTypeName(const std::string &colName, TTree *tree, RDataSource *ds, RDefineBase *define, bool vector2rvec)
Return a string containing the type of the given branch.
std::shared_ptr< RNodeBase > UpcastNode(std::shared_ptr< RNodeBase > ptr)
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...
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:
bool IsInternalColumn(std::string_view colName)
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)
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)
RInterface<::ROOT::Detail::RDF::RNodeBase, void > RNode
ROOT type_traits extensions.
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...
void EnableImplicitMT(UInt_t numthreads=0)
Enable ROOT's implicit multi-threading for all objects and methods that provide an internal paralleli...
Bool_t IsImplicitMTEnabled()
Returns true if the implicit multi-threading in ROOT is enabled.
void DisableImplicitMT()
Disables the implicit multi-threading in ROOT (see EnableImplicitMT).
A collection of options to steer the creation of the dataset on file.
bool fLazy
Do not start the event loop when Snapshot is called.
A struct which stores the parameters of a TH1D.
std::shared_ptr<::TH1D > GetHistogram() const
A struct which stores the parameters of a TH2D.
std::shared_ptr<::TH2D > GetHistogram() const
A struct which stores the parameters of a TH3D.
std::shared_ptr<::TH3D > GetHistogram() const
A struct which stores the parameters of a TProfile.
std::shared_ptr<::TProfile > GetProfile() const
A struct which stores the parameters of a TProfile2D.
std::shared_ptr<::TProfile2D > GetProfile() const
Lightweight storage for a collection of types.