Skip to content

Do analytical Barlow Beeston before likelihood evaluation #1137

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
anigamova merged 2 commits into from
Dec 4, 2025
Merged

Do analytical Barlow Beeston before likelihood evaluation #1137

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
9204f16
Cache pointers CMSHistSum and ErrorPropagator in CachingNLL
guitargeek Dec 1, 2025
1e3f6e3
Do analytical Barlow Beeston before likelihood evaluation
guitargeek Oct 29, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
7 changes: 3 additions & 4 deletions interface/CMSHistErrorPropagator.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,6 @@ class CMSHistErrorPropagator : public RooAbsReal {
std::vector<double> x2;
std::vector<double> res;
std::vector<double> gobs;
std::set<RooAbsArg*> dirty_prop;
std::vector<RooRealVar*> push_res;
};
public:
Expand Down Expand Up @@ -77,7 +76,9 @@ class CMSHistErrorPropagator : public RooAbsReal {

friend class CMSHistV<CMSHistErrorPropagator>;

protected:
void runBarlowBeeston() const;

protected:
RooRealProxy x_;
RooListProxy funcs_;
RooListProxy coeffs_;
Expand Down Expand Up @@ -109,8 +110,6 @@ class CMSHistErrorPropagator : public RooAbsReal {
void initialize() const;
void updateCache(int eval = 1) const;

void runBarlowBeeston() const;


private:
ClassDefOverride(CMSHistErrorPropagator,1)
Expand Down
8 changes: 3 additions & 5 deletions interface/CMSHistSum.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,6 @@ class CMSHistSum : public RooAbsReal {
std::vector<double> x2;
std::vector<double> res;
std::vector<double> gobs;
std::set<RooAbsArg*> dirty_prop;
std::vector<RooRealVar*> push_res;
};
public:
Expand Down Expand Up @@ -80,7 +79,9 @@ class CMSHistSum : public RooAbsReal {

void injectExternalMorph(int idx, CMSExternalMorph& morph);

protected:
void runBarlowBeeston() const;

protected:
RooRealProxy x_;

RooListProxy morphpars_;
Expand Down Expand Up @@ -143,9 +144,6 @@ class CMSHistSum : public RooAbsReal {
void updateCache() const;
inline double smoothStepFunc(double x, int const& ip) const;


void runBarlowBeeston() const;

void updateMorphs() const;


Expand Down
5 changes: 5 additions & 0 deletions interface/CachingNLL.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,8 @@
#include "SimpleConstraintGroup.h"

class RooMultiPdf;
class CMSHistSum;
class CMSHistErrorPropagator;

// Part zero: ArgSet checker
namespace cacheutils {
Expand Down Expand Up @@ -127,6 +129,7 @@ class CachingAddNLL : public RooAbsReal {
void updateZeroPoint() { clearZeroPoint(); setZeroPoint(); }
void propagateData();
void setAnalyticBarlowBeeston(bool flag);
void runAnalyticBarlowBeeston();
/// note: setIncludeZeroWeights(true) won't have effect unless you also re-call setData
virtual void setIncludeZeroWeights(bool includeZeroWeights) ;
RooSetProxy & params() { return params_; }
Expand All @@ -145,6 +148,8 @@ class CachingAddNLL : public RooAbsReal {
mutable std::vector<std::unique_ptr<RooAbsReal>> prods_;
mutable std::vector<RooAbsReal*> integrals_;
mutable std::vector<std::pair<const RooMultiPdf*,CachingPdfBase*> > multiPdfs_;
mutable std::vector<CMSHistSum const*> histSums_;
mutable std::vector<CMSHistErrorPropagator const*> histErrorPropagators_;
Comment on lines +151 to +152
Copy link
Copy Markdown
Contributor

@coderabbitai coderabbitai bot Oct 29, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ Potential issue | 🔴 Critical

🧩 Analysis chain

Non‑owning raw pointers: verify population and clearing.

Ensure these are:

  • Populated once during setup and deduplicated.
  • Cleared on data/PDF reset to avoid stale pointers.

Run:

🏁 Script executed:

#!/bin/bash
# Where are the vectors filled/cleared and the new method called?
rg -nP 'histSums_|histErrorPropagators_' -C2
rg -n 'runAnalyticBarlowBeeston\\(' -C3

Length of output: 2296

🏁 Script executed:

#!/bin/bash
# Search for clearing/resetting these vectors
rg -nP 'histSums_\s*\.|histErrorPropagators_\s*\.' -g '*.cc' -g '*.h' | grep -E '(clear|erase|reset|resize)'
# Look for methods that might reset cached data
rg -nP '(reset|clear|flush)' -A5 -B2 interface/CachingNLL.h src/CachingNLL.cc | head -100
# Search for where these vectors are modified/cleared in context
rg -n 'histSums_\s*=|histErrorPropagators_\s*=' src/CachingNLL.cc interface/CachingNLL.h

Length of output: 241

🏁 Script executed:

#!/bin/bash
# Find the files first
fd -e h -e cc 'CachingNLL' 

# If found, search for clearing operations
fd -e cc 'CachingNLL' -x rg -n 'clear|erase' {} \;
fd -e h 'CachingNLL' -x rg -n 'clear|erase' {} \;

# Search for destructor or reset methods
fd -e cc 'CachingNLL' -x rg -n '~CachingAddNLL|reset.*Data' {} \;

Length of output: 1293

🏁 Script executed:

#!/bin/bash
# Specifically search entire repo for histSums_ and histErrorPropagators_
rg -n 'histSums_' 
rg -n 'histErrorPropagators_' 

# Check the context around line 525 where integrals_ is cleared with other vectors
fd -e cc 'CachingNLL' -x rg -n -B5 -A5 'integrals_\.clear' {} \;

# Look for the method that clears integrals_ and other vectors
fd -e cc 'CachingNLL' -x rg -n -B10 'integrals_\.clear.*pdfs_\.clear' {} \;

Length of output: 1871

Confirm: Non-owning raw pointers must be cleared on reset.

The vectors histSums_ and histErrorPropagators_ are populated during PDF setup (src/CachingNLL.cc:588, 591) but are never cleared. This is inconsistent with other mutable cache vectors (integrals_, pdfs_, coeffs_, multiPdfs_) which are explicitly cleared in setup_() (line 525) and the destructor (line 413). Without clearing, these vectors retain stale pointers when PDFs or data are reset, creating a use-after-free risk in methods like propagateData(), setAnalyticBarlowBeeston(), and runAnalyticBarlowBeeston().

  • Add histSums_.clear() and histErrorPropagators_.clear() in setup_() after line 525 and in the destructor after line 413.
  • Verify deduplication logic during population (lines 587–592) to prevent duplicate pointers.
🤖 Prompt for AI Agents 
In interface/CachingNLL.h around lines 152–153, the mutable vectors histSums_
and histErrorPropagators_ hold non-owning raw pointers that are never cleared;
add histSums_.clear() and histErrorPropagators_.clear() in setup_() immediately
after the other cache clears (after line 525) and likewise in the destructor
after the other clears (after line 413) to avoid stale pointers, and while here,
inspect and fix the deduplication logic used when populating these vectors
(around src/CachingNLL.cc lines 587–592) to ensure duplicate pointers are not
inserted (e.g., check for existing pointer before push_back).

mutable std::vector<Double_t> partialSum_;
mutable std::vector<Double_t> workingArea_;
mutable bool isRooRealSum_, fastExit_;
Expand Down
12 changes: 1 addition & 11 deletions src/CMSHistErrorPropagator.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -99,7 +99,6 @@ void CMSHistErrorPropagator::initialize() const {
initialized_ = true;
}


void CMSHistErrorPropagator::updateCache(int eval) const {
initialize();

Expand Down Expand Up @@ -153,7 +152,6 @@ void CMSHistErrorPropagator::updateCache(int eval) const {


if (!binsentry_.good() || eval != last_eval_) {
runBarlowBeeston();
// bintypes might have size == 0 if we never ran setupBinPars()
for (unsigned j = 0; j < bintypes_.size(); ++j) {
cache_[j] = valsum_[j];
Expand Down Expand Up @@ -192,9 +190,7 @@ void CMSHistErrorPropagator::updateCache(int eval) const {
}

void CMSHistErrorPropagator::runBarlowBeeston() const {
if (!bb_.init) return;
RooAbsArg::setDirtyInhibit(true);

updateCache();
const unsigned n = bb_.use.size();
for (unsigned j = 0; j < n; ++j) {
bb_.dat[j] = data_[bb_.use[j]];
Expand All @@ -215,10 +211,6 @@ void CMSHistErrorPropagator::runBarlowBeeston() const {
for (unsigned j = 0; j < n; ++j) {
if (toterr_[bb_.use[j]] > 0.) bb_.push_res[j]->setVal(bb_.res[j]);
}
RooAbsArg::setDirtyInhibit(false);
for (RooAbsArg *arg : bb_.dirty_prop) {
arg->setValueDirty();
}
}

void CMSHistErrorPropagator::setAnalyticBarlowBeeston(bool flag) const {
Expand All @@ -240,7 +232,6 @@ void CMSHistErrorPropagator::setAnalyticBarlowBeeston(bool flag) const {
bb_.x2.clear();
bb_.res.clear();
bb_.gobs.clear();
bb_.dirty_prop.clear();
bb_.push_res.clear();
bb_.init = false;
}
Expand All @@ -254,7 +245,6 @@ void CMSHistErrorPropagator::setAnalyticBarlowBeeston(bool flag) const {
// std::cout << "Skipping " << this << " " << this->GetName() << "\n";
} else {
// std::cout << "Adding " << arg << " " << arg->GetName() << "\n";
bb_.dirty_prop.insert(arg);
auto as_gauss = dynamic_cast<RooGaussian*>(arg);
if (as_gauss) {
auto gobs = dynamic_cast<RooAbsReal*>(as_gauss->findServer(TString(vbinpars_[j][0]->GetName())+"_In"));
Expand Down
15 changes: 1 addition & 14 deletions src/CMSHistSum.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -336,7 +336,6 @@ inline double CMSHistSum::smoothStepFunc(double x, int const& ip) const {
return 0.125 * xnorm * (xnorm2 * (3. * xnorm2 - 10.) + 15);
}


void CMSHistSum::updateCache() const {
initialize();

Expand Down Expand Up @@ -379,10 +378,6 @@ void CMSHistSum::updateCache() const {


if (!binsentry_.good()) {
#if HFVERBOSE > 0
std::cout << "Calling runBarlowBeeston\n";
#endif
runBarlowBeeston();
// bintypes might have size == 0 if we never ran setupBinPars()
#if HFVERBOSE > 0
std::cout << "Assigning bin shifts\n";
Expand Down Expand Up @@ -418,9 +413,7 @@ void CMSHistSum::updateCache() const {
}

void CMSHistSum::runBarlowBeeston() const {
if (!bb_.init) return;
RooAbsArg::setDirtyInhibit(true);

updateCache();
const unsigned n = bb_.use.size();
for (unsigned j = 0; j < n; ++j) {
bb_.dat[j] = data_[bb_.use[j]];
Expand All @@ -441,10 +434,6 @@ void CMSHistSum::runBarlowBeeston() const {
for (unsigned j = 0; j < n; ++j) {
if (toterr_[bb_.use[j]] > 0.) bb_.push_res[j]->setVal(bb_.res[j]);
}
RooAbsArg::setDirtyInhibit(false);
for (RooAbsArg *arg : bb_.dirty_prop) {
arg->setValueDirty();
}
}

void CMSHistSum::setAnalyticBarlowBeeston(bool flag) const {
Expand All @@ -466,7 +455,6 @@ void CMSHistSum::setAnalyticBarlowBeeston(bool flag) const {
bb_.x2.clear();
bb_.res.clear();
bb_.gobs.clear();
bb_.dirty_prop.clear();
bb_.push_res.clear();
bb_.init = false;
}
Expand All @@ -480,7 +468,6 @@ void CMSHistSum::setAnalyticBarlowBeeston(bool flag) const {
// std::cout << "Skipping " << this << " " << this->GetName() << "\n";
} else {
// std::cout << "Adding " << arg << " " << arg->GetName() << "\n";
bb_.dirty_prop.insert(arg);
auto as_gauss = dynamic_cast<RooGaussian*>(arg);
if (as_gauss) {
auto gobs = dynamic_cast<RooAbsReal*>(as_gauss->findServer(TString(vbinpars_[j][0]->GetName())+"_In"));
Expand Down
66 changes: 43 additions & 23 deletions src/CachingNLL.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -568,13 +568,16 @@

multiPdfs_.clear();
for (auto &itp : pdfs_) {
bool isMultiPdf = itp->pdf()->IsA()->InheritsFrom(RooMultiPdf::Class());
if (isMultiPdf) {
const RooMultiPdf *mpdf = dynamic_cast<const RooMultiPdf*>(itp->pdf());
multiPdfs_.push_back(std::make_pair(mpdf, itp.get()));
}
if (auto* mpdf = dynamic_cast<const RooMultiPdf*>(itp->pdf())) {
multiPdfs_.emplace_back(mpdf, itp.get());
}
if (auto* hist = dynamic_cast<const CMSHistSum*>(itp->pdf())) {
histSums_.emplace_back(hist);
}
if (auto* hist = dynamic_cast<const CMSHistErrorPropagator*>(itp->pdf())) {
histErrorPropagators_.emplace_back(hist);
}
}
Comment on lines +571 to 580
Copy link
Copy Markdown
Contributor

@coderabbitai coderabbitai bot Dec 1, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ Potential issue | 🟠 Major

Clear vectors before populating to prevent stale/duplicate pointers.

When setup_() is called (e.g., after data changes), histSums_ and histErrorPropagators_ are populated without first clearing them, while other vectors like integrals_, pdfs_, coeffs_ are cleared at line 512. This can cause duplicate entries on re-setup or stale pointers if PDFs change.

Apply this diff:

 for (int i = 0, n = integrals_.size(); i < n; ++i) delete integrals_[i];
 integrals_.clear(); pdfs_.clear(); coeffs_.clear(); prods_.clear();
+histSums_.clear(); histErrorPropagators_.clear();

Committable suggestion skipped: line range outside the PR's diff.

🤖 Prompt for AI Agents 
In src/CachingNLL.cc around lines 571-580, histSums_ and histErrorPropagators_
are appended to without being cleared on re-setup, causing duplicate or stale
entries; before the loop that populates these vectors (near the existing clears
at line ~512), call histSums_.clear(); and histErrorPropagators_.clear(); (or
assign empty vectors) so those containers are emptied prior to repopulation.


}

void
Expand All @@ -592,6 +595,9 @@
#ifdef DEBUG_CACHE
PerfCounter::add("CachingAddNLL::evaluate called");
#endif
// The very first thing we do before any evaluation: run the analytical
// minimization of Barlow-Beeston nuisance parameters.
const_cast<CachingAddNLL&>(*this).runAnalyticBarlowBeeston();

std::fill( partialSum_.begin(), partialSum_.end(), 0.0 );

Expand Down Expand Up @@ -789,27 +795,30 @@
}

void cacheutils::CachingAddNLL::propagateData() {
for (auto const& funci : pdfs_) {
if ( auto pdf = dynamic_cast<CMSHistErrorPropagator const*>(funci->pdf()); pdf != nullptr ) {
pdf->setData(*data_);
}
else if ( auto pdf = dynamic_cast<CMSHistSum const*>(funci->pdf()); pdf != nullptr ) {
pdf->setData(*data_);
}
}
for (auto* hist : histErrorPropagators_) {
hist->setData(*data_);
}
for (auto* hist : histSums_) {
hist->setData(*data_);
}
}


void cacheutils::CachingAddNLL::setAnalyticBarlowBeeston(bool flag) {
for (auto const& funci : pdfs_) {
if ( auto pdf = dynamic_cast<CMSHistErrorPropagator const*>(funci->pdf()); pdf != nullptr ) {
pdf->setAnalyticBarlowBeeston(flag);
}
if ( auto pdf = dynamic_cast<CMSHistSum const*>(funci->pdf()); pdf != nullptr ) {
pdf->setAnalyticBarlowBeeston(flag);
}
for (auto* hist : histErrorPropagators_) {
hist->setAnalyticBarlowBeeston(flag);
}
for (auto* hist : histSums_) {
hist->setAnalyticBarlowBeeston(flag);
}
}

}
void cacheutils::CachingAddNLL::runAnalyticBarlowBeeston() {
for (auto* hist : histErrorPropagators_) {
hist->runBarlowBeeston();
}
for (auto* hist : histSums_) {
hist->runBarlowBeeston();
}
}

cacheutils::CachingSimNLL::CachingSimNLL(RooSimultaneous *pdf, RooAbsData *data, const RooArgSet *nuis) :
Expand Down Expand Up @@ -1010,6 +1019,17 @@
#ifdef DEBUG_CACHE
PerfCounter::add("CachingSimNLL::evaluate called");
#endif

// The very first thing we do before any evaluation: run the analytical
// minimization of Barlow-Beeston nuisance parameters.
for (size_t i = 0; i < pdfs_.size(); ++i) {
if (!pdfs_[i])
continue;

Check warning on line 1027 in src/CachingNLL.cc

View check run for this annotation

Codecov / codecov/patch

src/CachingNLL.cc#L1027 

Added line #L1027 was not covered by tests
if (!channelMasks_.empty() && channelMasks_[i]->getVal() != 0.)
continue;
pdfs_[i]->runAnalyticBarlowBeeston();
}

static bool gentleNegativePenalty_ = runtimedef::get("GENTLE_LEE");
DefaultAccumulator<double> ret = 0;
for (std::size_t idx = 0; idx < pdfs_.size(); ++idx) {
Expand Down
Loading