Merge pull request #164 from gpetruc/slc6-root5.34.17

merge updates from development branch to SL6

interface/CloseCoutSentry.h

@@ -13,13 +13,15 @@ class CloseCoutSentry {
         // break through any sentry, even the ones above myself (for critical error messages, or debug)
         static void breakFree() ;
         FILE *trueStdOut();
+        static FILE *trueStdOutGlobal();
     private:
         bool silent_;
         static int fdOut_, fdErr_;
         static bool open_;
         // always clear, even if I was not the one closing it
         void static reallyClear() ;
         static FILE *trueStdOut_; 
+        static CloseCoutSentry *owner_;
         bool stdOutIsMine_;
 };
 

interface/utils.h

@@ -3,6 +3,7 @@
 
 #include <vector>
 #include <string>
+#include <unordered_map>
 struct RooDataHist;
 struct RooAbsData;
 struct RooAbsPdf;
@@ -93,8 +94,8 @@ namespace utils {
     // Set range of physics model parameters
     void setModelParameterRanges( const std::string & setPhysicsModelParameterRangeExpression, const RooArgSet & params);
 
-    bool checkParameterBoundary( const RooRealVar &);
-    bool checkParameterBoundaries( const RooArgSet &);
+    bool isParameterAtBoundary( const RooRealVar &);
+    bool anyParameterAtBoundaries( const RooArgSet &, int verbosity);
 
     void reorderCombinations(std::vector<std::vector<int> > &, const std::vector<int> &, const std::vector<int> &);
     std::vector<std::vector<int> > generateCombinations(const std::vector<int> &vec);

src/CascadeMinimizer.cc

@@ -283,10 +283,11 @@ bool CascadeMinimizer::minimize(int verbose, bool cascade)
     nllParams->remove(CascadeMinimizerGlobalConfigs::O().pdfCategories);
     nllParams->remove(CascadeMinimizerGlobalConfigs::O().parametersOfInterest);
 
-    bool boundariesOk = utils::checkParameterBoundaries(*nllParams);
-    if(!boundariesOk){
-      std::cout << " [WARNING] After the fit some parameters are at their boundary (see above)." << std::endl;
-      std::cout << " [WARNING] Are you sure your model is correct?" << std::endl;
+    bool boundariesNotOk = utils::anyParameterAtBoundaries(*nllParams, verbose);
+    if(boundariesNotOk && verbose >= 1){
+      fprintf(CloseCoutSentry::trueStdOutGlobal(),
+        " [WARNING] After the fit some parameters are at their boundary.\n"
+        " [WARNING] Are you sure your model is correct?\n");
     }
 
     return ret;

src/CloseCoutSentry.cc

@@ -7,6 +7,8 @@ bool CloseCoutSentry::open_ = true;
 int  CloseCoutSentry::fdOut_ = 0;
 int  CloseCoutSentry::fdErr_ = 0;
 FILE * CloseCoutSentry::trueStdOut_ = 0;
+CloseCoutSentry *CloseCoutSentry::owner_ = 0;
+
 
 CloseCoutSentry::CloseCoutSentry(bool silent) :
     silent_(silent), stdOutIsMine_(false)
@@ -20,6 +22,7 @@ CloseCoutSentry::CloseCoutSentry(bool silent) :
             }
             freopen("/dev/null", "w", stdout);
             freopen("/dev/null", "w", stderr);
+            owner_ = this;
         } else {
             silent_ = false; 
         }
@@ -47,6 +50,7 @@ void CloseCoutSentry::reallyClear()
         sprintf(buf, "/dev/fd/%d", fdOut_); freopen(buf, "w", stdout);
         sprintf(buf, "/dev/fd/%d", fdErr_); freopen(buf, "w", stderr);
         open_   = true;
+        owner_ = 0;
         fdOut_ = fdErr_ = 0; 
     }
 }
@@ -56,6 +60,12 @@ void CloseCoutSentry::breakFree()
     reallyClear();
 }
 
+FILE *CloseCoutSentry::trueStdOutGlobal()
+{
+    if (!owner_) return stdout;
+    return owner_->trueStdOut();
+}
+
 FILE *CloseCoutSentry::trueStdOut() 
 {
     if (open_) return stdout;

src/utils.cc

@@ -6,6 +6,7 @@
 #include <sstream>
 #include <cmath>
 #include <vector>
+#include <unordered_map>
 #include <string>
 #include <memory>
 #include <typeinfo>
@@ -414,9 +415,11 @@ void utils::copyAttributes(const RooAbsArg &from, RooAbsArg &to) {
     if (!attribs.empty()) {
         for (std::set<std::string>::const_iterator it = attribs.begin(), ed = attribs.end(); it != ed; ++it) to.setAttribute(it->c_str());
     }
-    const std::map<std::string, std::string> strattribs = from.stringAttributes();
+    const std::map
+<std::string, std::string> strattribs = from.stringAttributes();
     if (!strattribs.empty()) {
-        for (std::map<std::string,std::string>::const_iterator it = strattribs.begin(), ed = strattribs.end(); it != ed; ++it) to.setStringAttribute(it->first.c_str(), it->second.c_str());
+        for (std::map
+        <std::string,std::string>::const_iterator it = strattribs.begin(), ed = strattribs.end(); it != ed; ++it) to.setStringAttribute(it->first.c_str(), it->second.c_str());
     }
 }
 
@@ -695,7 +698,8 @@ std::vector<std::vector<int> > utils::generateCombinations(const std::vector<int
 }
 
 
-bool utils::checkParameterBoundary( const RooRealVar &param ){
+bool utils::isParameterAtBoundary( const RooRealVar &param ){
+
     double vMin = param.getMin();
     double vMax = param.getMax();
     double val = param.getVal();
@@ -708,25 +712,40 @@ bool utils::checkParameterBoundary( const RooRealVar &param ){
     float nSigma=1.0;
 
     if(pullMin < nSigma || pullMax < nSigma){
-        CloseCoutSentry::breakFree();
-        std::cout << "  [WARNING] Found "<<param.GetName()<< " at " << std::min(pullMin,pullMax) << " sigma of one of its boundaries:" << std::endl;
-        std::cout << "   "; param.Print();
-        return false;
+        return true;
     }
-
-    return true;
+    
+    return false;
 }
 
-bool utils::checkParameterBoundaries( const RooArgSet &params ){
-
-    bool isNoneBad = true;
+
+bool utils::anyParameterAtBoundaries( const RooArgSet &params, int verbosity ){
+
+    static std::unordered_map<std::string, unsigned char> timesFoundAtBoundary;
+    bool isAnyBad = false;
 
     RooLinkedListIter iter = params.iterator(); int i = 0;
     for (RooRealVar *a = (RooRealVar *) iter.Next(); a != 0; a = (RooRealVar *) iter.Next(), ++i) {
-        bool isBad = checkParameterBoundary(*a);
-        isNoneBad &= isBad;
+
+        bool isBad = isParameterAtBoundary(*a);
+
+        if(isBad){
+            std::string varName((*a).GetName());
+
+            if( verbosity >= 9 || (timesFoundAtBoundary[varName] < 3 && verbosity > -1) ){
+                fprintf(CloseCoutSentry::trueStdOutGlobal(),"  [WARNING] Found [%s] at boundary. \n", (*a).GetName());
+                std::cout << "       "; (*a).Print();
+            }
+
+            timesFoundAtBoundary[varName]++;
+        }
+
+        isAnyBad |= isBad;
     }
 
-    return isNoneBad;
+    // for( std::unordered_map<std::string, unsigned char>::value_type e : timesFoundAtBoundary ){
+    //     printf("e %s -> %i\n", e.first.c_str(), e.second);
+    // }
+
+    return isAnyBad;
 }
-

test/diffNuisances.py

@@ -186,6 +186,18 @@
 elif options.format == "html":
     print "</table></body></html>"
 
+def getGraph(hist,shift):
+
+   gr = ROOT.TGraphErrors()
+   gr.SetName(hist.GetName())
+   for i in range(hist.GetNbinsX()):
+     x = hist.GetBinCenter(i+1)+shift
+     y = hist.GetBinContent(i+1)
+     e = hist.GetBinError(i+1)
+     gr.SetPoint(i,x,y)
+     gr.SetPointError(i,float(abs(shift))*0.8,e)
+   return gr
+
 if options.plotfile:
     import ROOT
     fout = ROOT.TFile(options.plotfile,"RECREATE")
@@ -205,6 +217,10 @@
     fout.WriteTObject(canvas)
 
     canvas_nuis = ROOT.TCanvas("nuisancs", "nuisances", 900, 600)
+    hist_fit_e_s = hist_fit_s.Clone()
+    hist_fit_e_b = hist_fit_b.Clone()
+    hist_fit_s = getGraph(hist_fit_s,-0.1)
+    hist_fit_b = getGraph(hist_fit_b, 0.1)
     hist_fit_s.SetLineColor(ROOT.kRed)
     hist_fit_s.SetMarkerColor(ROOT.kRed)
     hist_fit_b.SetLineColor(ROOT.kBlue)
@@ -219,11 +235,15 @@
     hist_prefit.SetTitle("Nuisance Paramaeters")
     hist_prefit.SetLineColor(ROOT.kBlack)
     hist_prefit.SetFillColor(ROOT.kGray)
+    hist_prefit.SetMaximum(3)
+    hist_prefit.SetMinimum(-3)
     hist_prefit.Draw("E2")
     hist_prefit.Draw("histsame")
-    hist_fit_b.Draw("E1Psame")
-    hist_fit_s.Draw("E1Psame")
+    hist_fit_b.Draw("EPsame")
+    hist_fit_s.Draw("EPsame")
+    canvas_nuis.SetGridx()
     canvas_nuis.RedrawAxis()
+    canvas_nuis.RedrawAxis('g')
     leg=ROOT.TLegend(0.6,0.7,0.89,0.89)
     leg.SetFillColor(0)
     leg.SetTextFont(42)