From 694e73711fb821db56038eaa1752dc1014ce1ed7 Mon Sep 17 00:00:00 2001 From: Nick Wardle Date: Wed, 17 Jul 2024 13:39:34 +0200 Subject: [PATCH] Adding getProcessNorms() to CMSHistErrorPropagator --- interface/CMSHistErrorPropagator.h | 2 ++ src/CMSHistErrorPropagator.cc | 35 ++++++++++++++++++++++++ test/printWorkspaceNormalisations.py | 40 ++++++++++++++++------------ 3 files changed, 60 insertions(+), 17 deletions(-) diff --git a/interface/CMSHistErrorPropagator.h b/interface/CMSHistErrorPropagator.h index 35e3196038c..c9e746a680e 100644 --- a/interface/CMSHistErrorPropagator.h +++ b/interface/CMSHistErrorPropagator.h @@ -72,6 +72,8 @@ class CMSHistErrorPropagator : public RooAbsReal { RooArgList wrapperList() const; RooArgList const& coefList() const { return coeffs_; } RooArgList const& funcList() const { return funcs_; } + + std::map getProcessNorms() const; friend class CMSHistV; diff --git a/src/CMSHistErrorPropagator.cc b/src/CMSHistErrorPropagator.cc index cf9f27ea8fa..1871c67268f 100644 --- a/src/CMSHistErrorPropagator.cc +++ b/src/CMSHistErrorPropagator.cc @@ -549,5 +549,40 @@ RooArgList CMSHistErrorPropagator::wrapperList() const { return result; } +std::map CMSHistErrorPropagator::getProcessNorms() const { + + std::map vals_; + RooArgList clist(coefList()); + RooArgList plist(funcList()); + /*if (plist.getSize() == 1) { + CMSHistErrorPropagator *err = dynamic_cast(plist.at(0)); + if (err) { + clist.removeAll(); + plist.removeAll(); + clist.add(err->coefList()); + plist.add(err->wrapperList()); + } + } + */ + for (int i = 0, n = clist.getSize(); i < n; ++i) { + RooAbsReal *coeff = (RooAbsReal *) clist.at(i); + std::string coeffName = coeff->GetName(); + RooAbsReal* shape = (RooAbsReal*)plist.at(i); + std::unique_ptr myobs(shape->getObservables(*x_)); + TString normProdName = TString::Format("%s", coeff->GetName()); + RooAbsReal * normProd = nullptr; + if (coeff->ownedComponents()) { + normProd = dynamic_cast(coeff->ownedComponents()->find(normProdName)); + } + if (!normProd) { + RooAbsReal* integral = shape->createIntegral(*myobs); + normProd = new RooProduct(normProdName, "", RooArgList(*integral, *coeff)); + normProd->addOwnedComponents(*integral); + coeff->addOwnedComponents(*normProd); + } + vals_[normProdName.Data()] = normProd->getVal(); + } + return vals_; +} #undef HFVERBOSE diff --git a/test/printWorkspaceNormalisations.py b/test/printWorkspaceNormalisations.py index 90269d4c773..c3f6431a7c3 100644 --- a/test/printWorkspaceNormalisations.py +++ b/test/printWorkspaceNormalisations.py @@ -49,14 +49,14 @@ type="float", help="Only print values if yield is less than this threshold.", ) -parser.add_option( - "", - "--use-cms-histsum", - dest="use_cms_histsum", - default=False, - action="store_true", - help="Set to true if workspace built with --use-cms-histsum option.", -) +#parser.add_option( +# "", +# "--use-cms-histsum", +# dest="use_cms_histsum", +# default=False, +# action="store_true", +# help="Set to true if workspace built with --use-cms-histsum option.", +#) parser.add_option( "", "--procFilter", @@ -169,17 +169,23 @@ def find_chan_proc(name): else: chan_procs[chan] = [[proc, norm, 0, 0]] +use_cms_histsum = False + +all_props = ws.allFunctions().selectByName("prop_bin*") +if all_props.getSize()>0: use_cms_histsum=True + # Look for cases where chan stored as CMSHistSum, set flag -if options.use_cms_histsum: +if use_cms_histsum: chan_CMSHistSum_norms = {} - all_props = ws.allFunctions().selectByName("prop_bin*") + #all_props = ws.allFunctions().selectByName("prop_bin*") for chan in chan_procs.keys(): prop_it = all_props.createIterator() for i in range(all_props.getSize()): prop = prop_it.Next() prop_name = prop.GetName() if chan == prop_name.split("_bin")[-1]: - if type(prop) == ROOT.CMSHistSum: + types = [ROOT.CMSHistSum, ROOT.CMSHistErrorPropagator] + if type(prop) in types : chan_CMSHistSum_norms[chan] = dict(prop.getProcessNorms()) @@ -243,7 +249,7 @@ def printEndExpand(): continue print("
") print("Top-level normalization for process {proc0} -> {proc1_name}
".format(proc0=proc[0], proc1_name=proc[1].GetName())) - if options.use_cms_histsum: + if use_cms_histsum: if chan in chan_CMSHistSum_norms: default_val = chan_CMSHistSum_norms[chan][proc[1].GetName()] else: @@ -253,7 +259,7 @@ def printEndExpand(): default_norms[chan][proc[1].GetName()] = default_val if options.printValueOnly: - print("default value = {default_val}
".format(default_val=default_val)) + print("default value = {default_val:.5f}
".format(default_val=default_val)) else: if proc[2]: proc_norm_var = ws.function("n_exp_bin%s_proc_%s" % (chan, proc[0])) @@ -289,7 +295,7 @@ def printEndExpand(): print(" ... is a constant (RooRealVar)") printEndExpand() - print(" default value = ", default_val, "
") + print(" default value = %.5f "%default_val, "
") print("") print( @@ -320,7 +326,7 @@ def printEndExpand(): print("---------------------------------------------------------------------------") print(" Top-level normalization for process %s -> %s" % (proc[0], proc[1].GetName())) print("---------------------------------------------------------------------------") - if options.use_cms_histsum: + if use_cms_histsum: if chan in chan_CMSHistSum_norms: default_val = chan_CMSHistSum_norms[chan][proc[1].GetName()] else: @@ -330,7 +336,7 @@ def printEndExpand(): default_norms[chan][proc[1].GetName()] = default_val if options.printValueOnly: - print(" default value = ", default_val) + print(" default value = %.5f "%default_val) # if options.printValueOnly: print " --xcp %s:%s "%(chan,proc[0]), else: if proc[2]: @@ -354,7 +360,7 @@ def printEndExpand(): proc[1].Print() print(" ... is a constant (RooRealVar)") print(" -------------------------------------------------------------------------") - print(" default value = ", default_val) + print(" default value = %.5f "%default_val) # Save norms to json file if options.output_json != "":