Skip to content

Commit

Permalink
[RF] Mode numeric integration logic for codegen into MathFuncs
Browse files Browse the repository at this point in the history
  • Loading branch information
guitargeek committed Oct 28, 2024
1 parent 98582c0 commit 56d5efd
Show file tree
Hide file tree
Showing 2 changed files with 20 additions and 26 deletions.
18 changes: 18 additions & 0 deletions roofit/roofitcore/inc/RooFit/Detail/MathFuncs.h
Original file line number Diff line number Diff line change
Expand Up @@ -728,6 +728,24 @@ inline double bernsteinIntegral(double xlo, double xhi, double xmin, double xmax
return norm * (xmax - xmin);
}

template <class Function>
inline double simpleNumericIntegral(double *params, double const *xlArr, Function func, double xmin, double xmax)
{
const int n = 1000; // number of sampling points
double d = xmax - xmin;
double eps = d / n;
double out = 0.;
double obs[1];
for (int i = 0; i < n; ++i) {
obs[0] = xmin + eps * i;
double tmpA = func(params, obs, xlArr);
obs[0] = xmin + eps * (i + 1);
double tmpB = func(params, obs, xlArr);
out += (tmpA + tmpB) * 0.5 * eps;
}
return out;
}

} // namespace MathFuncs

} // namespace Detail
Expand Down
28 changes: 2 additions & 26 deletions roofit/roofitcore/src/RooRealIntegral.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -1048,37 +1048,13 @@ void RooRealIntegral::translate(RooFit::CodegenContext &ctx) const

auto &intVar = static_cast<RooAbsRealLValue &>(*_intList[0]);

std::string obsName = ctx.getTmpVarName();
std::string oldIntVarResult = ctx.getResult(intVar);
ctx.addResult(&intVar, "obs[0]");

std::string funcName = ctx.buildFunction(*_function, {});

std::stringstream ss;

ss << "double " << obsName << "[1];\n";

std::string resName = RooFit::Detail::makeValidVarName(GetName()) + "Result";
ctx.addResult(this, resName);
ctx.addToGlobalScope("double " + resName + " = 0.0;\n");

// TODO: once Clad has support for higher-order functions (follow also the
// Clad issue #637), we could refactor this code into an actual function
// instead of hardcoding it here as a string.
ss << "{\n"
<< " const int n = 1000; // number of sampling points\n"
<< " double d = " << intVar.getMax(intRange()) << " - " << intVar.getMin(intRange()) << ";\n"
<< " double eps = d / n;\n"
<< " for (int i = 0; i < n; ++i) {\n"
<< " " << obsName << "[0] = " << intVar.getMin(intRange()) << " + eps * i;\n"
<< " double tmpA = " << funcName << "(params, " << obsName << ", xlArr);\n"
<< " " << obsName << "[0] = " << intVar.getMin(intRange()) << " + eps * (i + 1);\n"
<< " double tmpB = " << funcName << "(params, " << obsName << ", xlArr);\n"
<< " " << resName << " += (tmpA + tmpB) * 0.5 * eps;\n"
<< " }\n"
<< "}\n";

ctx.addToGlobalScope(ss.str());
ctx.addResult(this, ctx.buildCall("RooFit::Detail::MathFuncs::simpleNumericIntegral", "params", "xlArr", funcName,
intVar.getMin(intRange()), intVar.getMax(intRange())));

ctx.addResult(&intVar, oldIntVarResult);
}
Expand Down

0 comments on commit 56d5efd

Please sign in to comment.