Skip to content

Commit

Permalink
[RF] Make it clear when to use RooSimultaneous::generateSimGlobal
Browse files Browse the repository at this point in the history
More details in the comments in the code added in this commit.

To make it short: the user probably wants to use `generateSimGlobal`
when not generating the index category value (or taking it from a proto
dataset), and `RooAbsPdf::generate()` otherwise.

To avoid confusion, there are now new errors when one does the opposite,
which is certainly a user error.
  • Loading branch information
guitargeek committed Feb 5, 2025
1 parent 74dd985 commit 81fa486
Show file tree
Hide file tree
Showing 3 changed files with 48 additions and 45 deletions.
48 changes: 38 additions & 10 deletions roofit/roofitcore/src/RooSimultaneous.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -976,16 +976,24 @@ RooAbsGenContext* RooSimultaneous::genContext(const RooArgSet &vars, const RooDa
RooArgSet catsAmongAllVars;
allVars.selectCommon(flattenedCatList(), catsAmongAllVars);

// Not generating index cat: return context for pdf associated with present index state
// Not generating index cat: we better error out because it's not clear what
// the user expects here. Does the user want to generate according to the
// currently-selected pdf? Or does the user want to generate global
// observable values according to the union of all category pdfs?
// Print an error and tell the user what to do to explicitly.
if(catsAmongAllVars.empty()) {
auto* proxy = static_cast<RooRealProxy*>(_pdfProxyList.FindObject(_indexCat->getCurrentLabel()));
if (!proxy) {
coutE(InputArguments) << "RooSimultaneous::genContext(" << GetName()
<< ") ERROR: no PDF associated with current state ("
<< _indexCat.arg().GetName() << "=" << _indexCat.arg().getCurrentLabel() << ")" << std::endl ;
return nullptr;
}
return static_cast<RooAbsPdf*>(proxy->absArg())->genContext(vars,prototype,auxProto,verbose) ;
coutE(InputArguments) << "RooSimultaneous::generateSimGlobal(" << GetName()
<< ") asking to generate without the index category!\n"
<< "It's not clear what to do. you probably want to either:\n"
<< "\n"
<< " 1. Generate according to the currently-selected pdf.\n"
<< " Please do this explicitly with:\n"
<< " simpdf->getPdf(simpdf->indexCat().getCurrentLabel())->generate(vars, ...)\n"
<< "\n"
<< " 1. Generate global observable values according to the union of all component pdfs.\n"
<< " For this, please use simpdf->generateSimGlobal(vars, ...)\n"
<< std::endl;
return nullptr;
}

RooArgSet catsAmongProtoVars;
Expand Down Expand Up @@ -1032,10 +1040,30 @@ RooDataHist* RooSimultaneous::fillDataHist(RooDataHist *hist,


////////////////////////////////////////////////////////////////////////////////
/// Special generator interface for generation of 'global observables' -- for RooStats tools
/// Special generator interface for generation of 'global observables' -- for RooStats tools.
///
/// \note Why one can't just use RooAbsPdf::generate()? That's becaues when
/// using the regular generate() method, a specific component pdf is selected
/// for each generated entry according to the index category value. However,
/// global observable values are independent of the current index category,
/// which can best be illustrated with the case where a global observable
/// corresponds to a nuisance parameter that is relevant for multiple channels.
/// So the interpretation of what is an entry in the generated dataset is very
/// different, hence the separate function.

RooFit::OwningPtr<RooDataSet> RooSimultaneous::generateSimGlobal(const RooArgSet& whatVars, Int_t nEvents)
{
// Generating the index category together with the global observables doesn't make any sense.
RooArgSet catsAmongAllVars;
whatVars.selectCommon(flattenedCatList(), catsAmongAllVars);
if(!catsAmongAllVars.empty()) {
coutE(InputArguments) << "RooSimultaneous::generateSimGlobal(" << GetName()
<< ") asking to generate global obserables at the same time as the index category!\n"
<< "This doesn't make any sense: global observables are generally not related to a specific channel.\n"
<< std::endl;
return nullptr;
}

// Make set with clone of variables (placeholder for output)
RooArgSet globClone;
whatVars.snapshot(globClone);
Expand Down
23 changes: 5 additions & 18 deletions tutorials/roofit/roostats/TwoSidedFrequentistUpperLimitWithBands.C
Original file line number Diff line number Diff line change
Expand Up @@ -317,24 +317,11 @@ void TwoSidedFrequentistUpperLimitWithBands(const char *infile = "", const char
}

// generate global observables
// need to be careful for simpdf.
// In ROOT 5.28 there is a problem with generating global observables
// with a simultaneous PDF. In 5.29 there is a solution with
// RooSimultaneous::generateSimGlobal, but this may change to
// the standard generate interface in 5.30.

RooSimultaneous *simPdf = dynamic_cast<RooSimultaneous *>(mc->GetPdf());
if (!simPdf) {
std::unique_ptr<RooDataSet> one{mc->GetPdf()->generate(*mc->GetGlobalObservables(), 1)};
const RooArgSet *values = one->get();
std::unique_ptr<RooArgSet> allVars{mc->GetPdf()->getVariables()};
allVars->assign(*values);
} else {
std::unique_ptr<RooDataSet> one{simPdf->generateSimGlobal(*mc->GetGlobalObservables(), 1)};
const RooArgSet *values = one->get();
std::unique_ptr<RooArgSet> allVars{mc->GetPdf()->getVariables()};
allVars->assign(*values);
}

std::unique_ptr<RooDataSet> one{mc->GetPdf()->generateSimGlobal(*mc->GetGlobalObservables(), 1)};
const RooArgSet *values = one->get();
std::unique_ptr<RooArgSet> allVars{mc->GetPdf()->getVariables()};
allVars->assign(*values);

// get test stat at observed UL in observed data
firstPOI->setVal(observedUL);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -293,23 +293,11 @@
toyData = mc.GetPdf().generate(mc.GetObservables(), Extended=True)

# generate global observables
# need to be careful for simpdf.
# In ROOT 5.28 there is a problem with generating global observables
# with a simultaneous PDF. In 5.29 there is a solution with
# RooSimultaneous::generateSimGlobal, but this may change to
# the standard generate interface in 5.30.

simPdf = ROOT.RooSimultaneous(mc.GetPdf())
if not simPdf:
one = mc.GetPdf().generate(mc.GetGlobalObservables(), 1)
values = one.get()
allVars = mc.GetPdf().getVariables()
allVars.assign(values)
else:
one = simPdf.generateSimGlobal(mc.GetGlobalObservables(), 1)
values = one.get()
allVars = mc.GetPdf().getVariables()
allVars.assign(values)

one = mc.GetPdf().generateSimGlobal(mc.GetGlobalObservables(), 1)
values = one.get()
allVars = mc.GetPdf().getVariables()
allVars.assign(values)

# get test stat at observed UL in observed data
firstPOI.setVal(observedUL)
Expand Down

0 comments on commit 81fa486

Please sign in to comment.