Skip to content

Commit a8ef8b0

Browse files
authored
[RF] Fix plotOn normalization for extended PDFs
This PR fixes an issue where the `plotOn` function incorrectly normalizes an extended PDF to the number of events in a dataset, rather than to its own predicted number of events. This behavior is misleading when the PDF's normalization is a key part of the model's prediction. The default logic in `RooAbsPdf::plotOn` is modified to check if a PDF is extended by calling `expectedEvents()`. If it is, the function now uses the PDF's prediction for normalization; otherwise, it retains the old behavior of scaling to the data for non-extended PDFs. Closes #18929.
1 parent a10d2af commit a8ef8b0

File tree

1 file changed

+15
-1
lines changed

1 file changed

+15
-1
lines changed

roofit/roofitcore/src/RooAbsPdf.cxx

Lines changed: 15 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -2145,7 +2145,21 @@ RooPlot* RooAbsPdf::plotOn(RooPlot* frame, RooLinkedList& cmdList) const
21452145
scaleFactor *= rangeNevt/nExpected ;
21462146

21472147
} else {
2148-
scaleFactor *= frame->getFitRangeNEvt()/nExpected ;
2148+
// First, check if the PDF *can* be extended.
2149+
if (this->canBeExtended()) {
2150+
// If it can, get the expected events.
2151+
const double nExp = expectedEvents(frame->getNormVars());
2152+
if (nExp > 0) {
2153+
// If the prediction is valid, use it for normalization.
2154+
scaleFactor *= nExp / nExpected;
2155+
} else {
2156+
// If prediction is not valid (e.g. 0), fall back to data.
2157+
scaleFactor *= frame->getFitRangeNEvt() / nExpected;
2158+
}
2159+
} else {
2160+
// If the PDF can't be extended, just use the data.
2161+
scaleFactor *= frame->getFitRangeNEvt() / nExpected;
2162+
}
21492163
}
21502164
} else if (stype==RelativeExpected) {
21512165
scaleFactor *= nExpected ;

0 commit comments

Comments
 (0)