From 24f6fad316b01a0e07d5df2f9c9fe27420272e8b Mon Sep 17 00:00:00 2001 From: Jonas Rembser Date: Wed, 19 Aug 2026 10:01:16 +0000 Subject: [PATCH] [RF] Fix counting Asimov generation with multiple floating parameters MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit When generating a counting Asimov dataset, the AsymptoticCalculator identified the expected value of an observable in a Gaussian or Poisson term as "the single non-constant server that is not the observable". This heuristic failed whenever both the mean and the width of a Gaussian were floating, could not see constness through derived quantities like RooFormulaVars, and silently set the observable to the value of the *sigma* parameter when the mean was constant but sigma floated. The heuristic is unnecessary: RooGaussian and RooPoisson expose their proxies via getX() and getMean(). Use the exact (x, mean) pair and set whichever of the two is the observable to the value of the other one. The width and the constness of the parameters no longer matter, and the direction where the mean is the observable (constraint terms with global observables) keeps working. The RooMultiVarGaussian path goes through the same helper with xVec()[i] and muVec()[i]. Fixes the four failure cases from ROOT-10096, covered by a new test. Fixes [ROOT-10096](https://its.cern.ch/jira/browse/ROOT-10069) 🤖 Done with the help of AI --- roofit/roostats/src/AsymptoticCalculator.cxx | 79 ++++++++----------- .../test/testAsymptoticCalculator.cxx | 38 +++++++++ 2 files changed, 72 insertions(+), 45 deletions(-) diff --git a/roofit/roostats/src/AsymptoticCalculator.cxx b/roofit/roostats/src/AsymptoticCalculator.cxx index 72cb592087252..b6d4e2ad935b4 100644 --- a/roofit/roostats/src/AsymptoticCalculator.cxx +++ b/roofit/roostats/src/AsymptoticCalculator.cxx @@ -69,8 +69,6 @@ The calculator can generate Asimov datasets from two kinds of PDFs: #include "TStopwatch.h" -#include - using namespace RooStats; using std::string, std::unique_ptr; @@ -939,41 +937,31 @@ void FillBins(const RooAbsPdf & pdf, const RooArgList &obs, RooAbsData & data, i } -bool setObsToExpected(std::span servers, const RooArgSet &obs, std::string const &errPrefix) +bool setObsToExpected(RooAbsArg &x, RooAbsArg &mean, const RooArgSet &obs, std::string const &errPrefix) { - RooRealVar *myobs = nullptr; - RooAbsReal *myexp = nullptr; - for (RooAbsArg *a : servers) { - if (obs.contains(*a)) { - if (myobs != nullptr) { - oocoutF(nullptr,Generation) << errPrefix << "Has two observables ?? " << std::endl; - return false; - } - myobs = dynamic_cast(a); - if (myobs == nullptr) { - oocoutF(nullptr,Generation) << errPrefix << "Observable is not a RooRealVar??" << std::endl; - return false; - } - } else { - if (!a->isConstant() ) { - if (myexp != nullptr) { - oocoutE(nullptr,Generation) << errPrefix << "Has two non-const arguments " << std::endl; - return false; - } - myexp = dynamic_cast(a); - if (myexp == nullptr) { - oocoutF(nullptr,Generation) << errPrefix << "Expected is not a RooAbsReal??" << std::endl; - return false; - } - } - } + // Figure out which of the two arguments is the observable that should be + // set to the expected value given by the other one. Usually the observable + // is "x", but also the mean parameter can be the observable: this happens + // for example in constraint terms, where the global observable takes the + // role of the mean. + const bool xIsObs = obs.contains(x); + const bool meanIsObs = obs.contains(mean); + if (xIsObs && meanIsObs) { + oocoutF(nullptr, Generation) << errPrefix << "Has two observables ?? " << std::endl; + return false; + } + if (!xIsObs && !meanIsObs) { + oocoutF(nullptr, Generation) << errPrefix << "No observable?" << std::endl; + return false; } - if (myobs == nullptr) { - oocoutF(nullptr,Generation) << errPrefix << "No observable?" << std::endl; + auto *myobs = dynamic_cast(xIsObs ? &x : &mean); + auto *myexp = dynamic_cast(xIsObs ? &mean : &x); + if (myobs == nullptr) { + oocoutF(nullptr, Generation) << errPrefix << "Observable is not a RooRealVar??" << std::endl; return false; } if (myexp == nullptr) { - oocoutF(nullptr,Generation) << errPrefix << "No observable?" << std::endl; + oocoutF(nullptr, Generation) << errPrefix << "Expected is not a RooAbsReal??" << std::endl; return false; } @@ -989,33 +977,34 @@ bool setObsToExpected(std::span servers, const RooArgSet &obs, std: //////////////////////////////////////////////////////////////////////////////// /// set observed value to the expected one /// works for Gaussian, Poisson or LogNormal -/// assumes mean parameter value is the argument not constant and not depending on observables -/// (if more than two arguments are not constant will use first one but print a warning !) /// need to iterate on the components of the Poisson to get n and nu (nu can be a RooAbsReal) /// (code from G. Petrucciani and extended by L.M.) -bool SetObsToExpected(RooAbsPdf &pdf, const RooArgSet &obs) +bool SetObsToExpected(RooGaussian &pdf, const RooArgSet &obs) { std::string const &errPrefix = "AsymptoticCalculator::SetObsExpected( " + std::string{pdf.ClassName()} + " ) : "; - std::vector servers; - for (RooAbsArg *a : pdf.servers()) { - servers.emplace_back(a); - } - return setObsToExpected(servers, obs, errPrefix); + return setObsToExpected(const_cast(pdf.getX()), const_cast(pdf.getMean()), obs, + errPrefix); +} + +bool SetObsToExpected(RooPoisson &pdf, const RooArgSet &obs) +{ + std::string const &errPrefix = "AsymptoticCalculator::SetObsExpected( " + std::string{pdf.ClassName()} + " ) : "; + return setObsToExpected(const_cast(pdf.getX()), const_cast(pdf.getMean()), obs, + errPrefix); } bool setObsToExpectedMultiVarGauss(RooMultiVarGaussian &mvgauss, const RooArgSet &obs) { // In the case of the multi-variate Gaussian, we need to iterate over the - // dimensions and treat the servers for each dimension separately. + // dimensions and treat the observable and mean for each dimension + // separately. std::string const &errPrefix = "AsymptoticCalculator::SetObsExpected( " + std::string{mvgauss.ClassName()} + " ) : "; - std::vector servers{nullptr, nullptr}; bool ret = true; for (std::size_t iDim = 0; iDim < mvgauss.xVec().size(); ++iDim) { - servers[0] = &mvgauss.xVec()[iDim]; - servers[1] = &mvgauss.muVec()[iDim]; - ret &= setObsToExpected(servers, obs, errPrefix + " : dim " + std::to_string(iDim) + " "); + ret &= setObsToExpected(mvgauss.xVec()[iDim], mvgauss.muVec()[iDim], obs, + errPrefix + " : dim " + std::to_string(iDim) + " "); } return ret; } diff --git a/roofit/roostats/test/testAsymptoticCalculator.cxx b/roofit/roostats/test/testAsymptoticCalculator.cxx index f49d8548f125a..da5defd4f19d1 100644 --- a/roofit/roostats/test/testAsymptoticCalculator.cxx +++ b/roofit/roostats/test/testAsymptoticCalculator.cxx @@ -109,3 +109,41 @@ TEST(AsymptoticCalculator, SignedTestStatistic) std::unique_ptr resRightSideCapped{calc.GetHypoTest()}; EXPECT_NEAR(resRightSideCapped->NullPValue(), resRightSide->NullPValue(), 1e-6); } + +// Check that counting Asimov datasets can be generated no matter which +// parameters are floating and even if the mean or width of a Gaussian are +// derived quantities (covers JIRA ROOT-10096). +TEST(AsymptoticCalculator, CountingAsimovDataSetFloatingParams) +{ + RooWorkspace ws; + ws.factory("obs[10.0, 0.0, 1000.0]"); + ws.factory("Poisson::poisson(obs, mean[20.0, 0.0, 1000.0])"); + ws.factory("Gaussian::gauss1(obs, mean, sigma[3.0, 1.0, 10.0])"); + ws.factory("expr::sqrt_mean('sqrt(@0)', mean)"); + ws.factory("Gaussian::gauss2(obs, mean, sqrt_mean)"); + ws.factory("expr::mean2('2 * @0', mean)"); + ws.factory("expr::sqrt_mean2('sqrt(@0)', mean2)"); + ws.factory("Gaussian::gauss3(obs, mean2, sqrt_mean2)"); + + RooArgSet observables{*ws.var("obs")}; + + auto checkAsimov = [&](const char *pdfName, double expectedObsVal) { + std::unique_ptr data{ + RooStats::AsymptoticCalculator::GenerateAsimovData(*ws.pdf(pdfName), observables)}; + ASSERT_NE(data, nullptr) << pdfName; + ASSERT_EQ(data->numEntries(), 1) << pdfName; + EXPECT_DOUBLE_EQ(data->get(0)->getRealValue("obs"), expectedObsVal) << pdfName; + }; + + checkAsimov("poisson", 20.0); + // Both mean and sigma floating: used to fail with "Has two non-const arguments". + checkAsimov("gauss1", 20.0); + // Width derived from the mean: also used to fail, with no workaround for gauss3. + checkAsimov("gauss2", 20.0); + checkAsimov("gauss3", 40.0); + + // With a constant mean and a floating sigma, the old server-based heuristic + // silently set the observable to the value of the sigma parameter. + ws.var("mean")->setConstant(true); + checkAsimov("gauss1", 20.0); +}