Skip to content
Snippets Groups Projects
Commit f2a41f41 authored by Jonas Rembser's avatar Jonas Rembser
Browse files

[RF] Correctly normalize max possible val in generation with proto data

The code branch without a proto data in `RooGenContext` already did
this, and now it's also done when using proto dataset generation.

Closes #12286.
parent 8c7de001
No related branches found
No related tags found
No related merge requests found
// Tests for the RooKeysPdf and friends
// Authors: Jonas Rembser, CERN 07/2022
#include <RooRealVar.h>
#include <RooGenericPdf.h>
#include <RooKeysPdf.h>
#include <RooNDKeysPdf.h>
#include <RooPlot.h>
#include <RooRealVar.h>
#include "gtest/gtest.h"
......@@ -54,3 +56,29 @@ TEST(RooKeysPdf, WeightedDataset)
xVal += 0.1;
}
}
// Test generation with proto data, covering GitHub issue #12286.
TEST(RooKeysPdf, GenerationWithProtoData)
{
using namespace RooFit;
RooRealVar x{"x", "", 0, 1};
RooGenericPdf pdfX{"pdf_x", "x", {x}};
std::unique_ptr<RooDataSet> dtBase{pdfX.generate(x, 10000)};
RooKeysPdf pdfKeys{"pdf_keys", "", x, *dtBase, RooKeysPdf::MirrorBoth};
RooRealVar y{"y", "", 0, 1};
RooDataSet proto{"proto_y", "", y};
proto.add(y);
std::unique_ptr<RooDataSet> dtKeysWithProto{pdfKeys.generate(x, NumEvents(10000), ProtoData(proto))};
std::unique_ptr<RooPlot> frame{x.frame()};
dtKeysWithProto->plotOn(frame.get());
pdfKeys.plotOn(frame.get());
// If the dataset generation worked, the chi-square is not too terrible
EXPECT_LE(frame->chiSquare(), 2.0);
}
......@@ -240,7 +240,8 @@ RooGenContext::RooGenContext(const RooAbsPdf &model, const RooArgSet &vars,
coutI(Generation) << "RooGenContext::ctor() prototype data provided, and "
<< "model supports analytical maximum finding in the full phase space: "
<< "can provide analytical pdf maximum to numeric generator" << endl ;
_maxVar = new RooRealVar("funcMax","function maximum",_pdfClone->maxVal(maxFindCode)) ;
double maxVal = _pdfClone->maxVal(maxFindCode) / _pdfClone->getNorm(&allVars);
_maxVar = new RooRealVar("funcMax", "function maximum", maxVal);
} else {
maxFindCode = _pdfClone->getMaxVal(_otherVars) ;
if (maxFindCode != 0) {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment