RooStats/HistFactory/Systematics.h

0001 // @(#)root/roostats:$Id$
0002 // Author: George Lewis, Kyle Cranmer
0003 /*************************************************************************
0004  * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers.               *
0005  * All rights reserved.                                                  *
0006  *                                                                       *
0007  * For the licensing terms see $ROOTSYS/LICENSE.                         *
0008  * For the list of contributors see $ROOTSYS/README/CREDITS.             *
0009  *************************************************************************/
0010
0011 #ifndef HISTFACTORY_SYSTEMATICS_H
0012 #define HISTFACTORY_SYSTEMATICS_H
0013
0014 #include <string>
0015 #include <fstream>
0016 #include <iostream>
0017
0018 #include "TH1.h"
0019 #include "RooStats/HistFactory/HistRef.h"
0020
0021 namespace RooStats{
0022 namespace HistFactory {
0023
0024   namespace Constraint {
0025     enum Type{ Gaussian, Poisson };
0026     std::string Name( Type type );
0027     Type GetType( const std::string& Name );
0028   }
0029
0030 /** \class OverallSys
0031  * \ingroup HistFactory
0032  * Configuration for a constrained overall systematic to scale sample normalisations.
0033  */
0034   class OverallSys {
0035
0036   public:
0037
0038     OverallSys() : fLow(0), fHigh(0) {}
0039
0040     void SetName( const std::string& Name ) { fName = Name; }
0041     const std::string& GetName() const { return fName; }
0042
0043     void SetLow( double Low )   { fLow  = Low; }
0044     void SetHigh( double High ) { fHigh = High; }
0045     double GetLow() const { return fLow; }
0046     double GetHigh() const { return fHigh; }
0047
0048     void Print(std::ostream& = std::cout) const;
0049     void PrintXML(std::ostream&) const;
0050
0051   protected:
0052     std::string fName;
0053     double fLow;
0054     double fHigh;
0055
0056   };
0057
0058 /** \class NormFactor
0059  * \ingroup HistFactory
0060  * Configuration for an \a un- constrained overall systematic to scale sample normalisations.
0061  */
0062   class NormFactor {
0063
0064   public:
0065
0066     NormFactor();
0067
0068     void SetName( const std::string& Name ) { fName = Name; }
0069     std::string GetName() const { return fName; }
0070
0071     void SetVal( double Val ) { fVal = Val; }
0072     double GetVal() const { return fVal; }
0073
0074     void SetLow( double Low )   { fLow  = Low; }
0075     void SetHigh( double High ) { fHigh = High; }
0076     double GetLow() const { return fLow; }
0077     double GetHigh() const { return fHigh; }
0078
0079     void Print(std::ostream& = std::cout) const;
0080     void PrintXML(std::ostream&) const;
0081
0082   protected:
0083
0084     std::string fName;
0085     double fVal;
0086     double fLow;
0087     double fHigh;
0088
0089   };
0090
0091
0092   /** ////////////////////////////////////////////////////////////////////////////////////////////
0093    * \class HistogramUncertaintyBase
0094    * \ingroup HistFactory
0095    * Base class to store the up and down variations for histogram uncertainties.
0096    * Use the derived classes for actual models.
0097    */
0098   class HistogramUncertaintyBase {
0099
0100   public:
0101
0102     HistogramUncertaintyBase() : fhLow(nullptr), fhHigh(nullptr) {}
0103     HistogramUncertaintyBase(const std::string& Name) : fName(Name), fhLow(nullptr), fhHigh(nullptr) {}
0104     HistogramUncertaintyBase(const HistogramUncertaintyBase& oth) :
0105     fName{oth.fName},
0106     fInputFileLow{oth.fInputFileLow}, fHistoNameLow{oth.fHistoNameLow}, fHistoPathLow{oth.fHistoPathLow},
0107     fInputFileHigh{oth.fInputFileHigh}, fHistoNameHigh{oth.fHistoNameHigh}, fHistoPathHigh{oth.fHistoPathHigh},
0108     fhLow{oth.fhLow ? static_cast<TH1*>(oth.fhLow->Clone()) : nullptr},
0109     fhHigh{oth.fhHigh ? static_cast<TH1*>(oth.fhHigh->Clone()) : nullptr} {
0110
0111     }
0112     HistogramUncertaintyBase(HistogramUncertaintyBase&&) = default;
0113
0114     virtual ~HistogramUncertaintyBase() {};
0115
0116
0117     // Need deep copies because the class owns its histograms.
0118     HistogramUncertaintyBase& operator=(const HistogramUncertaintyBase& oth) {
0119       fName = oth.fName;
0120       fInputFileLow = oth.fInputFileLow;
0121       fHistoNameLow = oth.fHistoNameLow;
0122       fHistoPathLow = oth.fHistoPathLow;
0123       fInputFileHigh = oth.fInputFileHigh;
0124       fHistoNameHigh = oth.fHistoNameHigh;
0125       fHistoPathHigh = oth.fHistoPathHigh;
0126       fhLow.reset(oth.fhLow ? static_cast<TH1*>(oth.fhLow->Clone()) : nullptr);
0127       fhHigh.reset(oth.fhHigh ? static_cast<TH1*>(oth.fhHigh->Clone()) : nullptr);
0128
0129       return *this;
0130     }
0131     HistogramUncertaintyBase& operator=(HistogramUncertaintyBase&&) = default;
0132
0133     virtual void Print(std::ostream& = std::cout) const;
0134     virtual void PrintXML(std::ostream&) const = 0;
0135     virtual void writeToFile( const std::string& FileName, const std::string& DirName );
0136
0137     void SetHistoLow(TH1* Low ) {Low->SetDirectory(nullptr); fhLow.reset(Low);}
0138     void SetHistoHigh(TH1* High ) {High->SetDirectory(nullptr); fhHigh.reset(High);}
0139
0140     const TH1* GetHistoLow() const {return fhLow.get();}
0141     const TH1* GetHistoHigh() const {return fhHigh.get();}
0142
0143     void SetName( const std::string& Name ) { fName = Name; }
0144     const std::string& GetName() const { return fName; }
0145
0146     void SetInputFileLow( const std::string& InputFileLow ) { fInputFileLow = InputFileLow; }
0147     void SetInputFileHigh( const std::string& InputFileHigh ) { fInputFileHigh = InputFileHigh; }
0148
0149     const std::string& GetInputFileLow() const { return fInputFileLow; }
0150     const std::string& GetInputFileHigh() const { return fInputFileHigh; }
0151
0152     void SetHistoNameLow( const std::string& HistoNameLow ) { fHistoNameLow = HistoNameLow; }
0153     void SetHistoNameHigh( const std::string& HistoNameHigh ) { fHistoNameHigh = HistoNameHigh; }
0154
0155     const std::string& GetHistoNameLow() const { return fHistoNameLow; }
0156     const std::string& GetHistoNameHigh() const { return fHistoNameHigh; }
0157
0158     void SetHistoPathLow( const std::string& HistoPathLow ) { fHistoPathLow = HistoPathLow; }
0159     void SetHistoPathHigh( const std::string& HistoPathHigh ) { fHistoPathHigh = HistoPathHigh; }
0160
0161     const std::string& GetHistoPathLow() const { return fHistoPathLow; }
0162     const std::string& GetHistoPathHigh() const { return fHistoPathHigh; }
0163
0164   protected:
0165
0166     std::string fName;
0167
0168     std::string fInputFileLow;
0169     std::string fHistoNameLow;
0170     std::string fHistoPathLow;
0171
0172     std::string fInputFileHigh;
0173     std::string fHistoNameHigh;
0174     std::string fHistoPathHigh;
0175
0176     // The Low and High Histograms
0177     std::unique_ptr<TH1> fhLow;
0178     std::unique_ptr<TH1> fhHigh;
0179
0180   };
0181
0182 /** \class HistoSys
0183  * \ingroup HistFactory
0184  * Configuration for a constrained, coherent shape variation of affected samples.
0185  */
0186 class HistoSys final : public HistogramUncertaintyBase {
0187 public:
0188   void PrintXML(std::ostream&) const override;
0189 };
0190
0191 /** \class HistoFactor
0192  * \ingroup HistFactory
0193  * Configuration for an *un*constrained, coherent shape variation of affected samples.
0194  */
0195   class HistoFactor final : public HistogramUncertaintyBase {
0196   public:
0197     void PrintXML(std::ostream&) const override;
0198   };
0199
0200 /** \class ShapeSys
0201  * \ingroup HistFactory
0202  * Constrained bin-by-bin variation of affected histogram.
0203  */
0204   class ShapeSys final : public HistogramUncertaintyBase {
0205
0206   public:
0207     ShapeSys() : fConstraintType(Constraint::Gaussian) {}
0208     ShapeSys(const ShapeSys& other) :
0209       HistogramUncertaintyBase(other),
0210       fConstraintType(other.fConstraintType) {}
0211     ShapeSys& operator=(const ShapeSys& oth) {
0212        if (this == &oth) return *this;
0213        HistogramUncertaintyBase::operator=(oth);
0214        fConstraintType = oth.fConstraintType;
0215        return *this;
0216     }
0217     ShapeSys& operator=(ShapeSys&&) = default;
0218
0219     void SetInputFile( const std::string& InputFile ) { fInputFileHigh = InputFile; }
0220     std::string GetInputFile() const { return fInputFileHigh; }
0221
0222     void SetHistoName( const std::string& HistoName ) { fHistoNameHigh = HistoName; }
0223     std::string GetHistoName() const { return fHistoNameHigh; }
0224
0225     void SetHistoPath( const std::string& HistoPath ) { fHistoPathHigh = HistoPath; }
0226     std::string GetHistoPath() const { return fHistoPathHigh; }
0227
0228     void Print(std::ostream& = std::cout) const override;
0229     void PrintXML(std::ostream&) const override;
0230     void writeToFile( const std::string& FileName, const std::string& DirName ) override;
0231
0232     const TH1* GetErrorHist() const {
0233       return fhHigh.get();
0234     }
0235     void SetErrorHist(TH1* hError) {
0236       fhHigh.reset(hError);
0237     }
0238
0239     void SetConstraintType( Constraint::Type ConstrType ) { fConstraintType = ConstrType; }
0240     Constraint::Type GetConstraintType() const { return fConstraintType; }
0241
0242   protected:
0243     Constraint::Type fConstraintType;
0244   };
0245
0246 /** \class ShapeFactor
0247  * \ingroup HistFactory
0248  * *Un*constrained bin-by-bin variation of affected histogram.
0249  */
0250   class ShapeFactor : public HistogramUncertaintyBase {
0251
0252   public:
0253     void Print(std::ostream& = std::cout) const override;
0254     void PrintXML(std::ostream&) const override;
0255     void writeToFile( const std::string& FileName, const std::string& DirName) override;
0256
0257     void SetInitialShape(TH1* shape) {
0258       fhHigh.reset(shape);
0259     }
0260     const TH1* GetInitialShape() const { return fhHigh.get(); }
0261
0262     void SetConstant(bool constant) { fConstant = constant; }
0263     bool IsConstant() const { return fConstant; }
0264
0265     bool HasInitialShape() const { return fHasInitialShape; }
0266
0267     void SetInputFile( const std::string& InputFile ) {
0268       fInputFileHigh = InputFile;
0269       fHasInitialShape=true;
0270     }
0271     const std::string& GetInputFile() const { return fInputFileHigh; }
0272
0273     void SetHistoName( const std::string& HistoName ) {
0274       fHistoNameHigh = HistoName;
0275       fHasInitialShape=true;
0276     }
0277     const std::string& GetHistoName() const { return fHistoNameHigh; }
0278
0279     void SetHistoPath( const std::string& HistoPath ) {
0280       fHistoPathHigh = HistoPath;
0281       fHasInitialShape=true;
0282     }
0283     const std::string& GetHistoPath() const { return fHistoPathHigh; }
0284
0285   protected:
0286
0287     bool fConstant = false;
0288
0289     // A histogram representing
0290     // the initial shape
0291     bool fHasInitialShape = false;
0292   };
0293
0294 /** \class StatError
0295  * \ingroup HistFactory
0296  * Statistical error of Monte Carlo predictions.
0297  */
0298   class StatError : public HistogramUncertaintyBase {
0299
0300   public:
0301     void Print(std::ostream& = std::cout) const override;
0302     void PrintXML(std::ostream&) const override;
0303     void writeToFile( const std::string& FileName, const std::string& DirName ) override;
0304
0305     void Activate( bool IsActive=true ) { fActivate = IsActive; }
0306     bool GetActivate() const { return fActivate; }
0307
0308     void SetUseHisto( bool UseHisto=true ) { fUseHisto = UseHisto; }
0309     bool GetUseHisto() const { return fUseHisto; }
0310
0311     void SetInputFile( const std::string& InputFile ) { fInputFileHigh = InputFile; }
0312     const std::string& GetInputFile() const { return fInputFileHigh; }
0313
0314     void SetHistoName( const std::string& HistoName ) { fHistoNameHigh = HistoName; }
0315     const std::string& GetHistoName() const { return fHistoNameHigh; }
0316
0317     void SetHistoPath( const std::string& HistoPath ) { fHistoPathHigh = HistoPath; }
0318     const std::string& GetHistoPath() const { return fHistoPathHigh; }
0319
0320
0321     const TH1* GetErrorHist() const {
0322       return fhHigh.get();
0323     }
0324     void SetErrorHist(TH1* Error) {
0325       fhHigh.reset(Error);
0326     }
0327
0328   protected:
0329
0330     bool fActivate = false;
0331     bool fUseHisto = false; // Use an external histogram for the errors
0332   };
0333
0334 /** \class StatErrorConfig
0335  * \ingroup HistFactory
0336  * Configuration to automatically assign nuisance parameters for the statistical
0337  * error of the Monte Carlo simulations.
0338  * The default is to assign a Poisson uncertainty to a bin when its statistical uncertainty
0339  * is larger than 5% of the bin content.
0340  */
0341   class StatErrorConfig {
0342
0343   public:
0344
0345     StatErrorConfig() : fRelErrorThreshold( .05 ), fConstraintType( Constraint::Poisson ) {;}
0346     void Print(std::ostream& = std::cout) const;
0347     void PrintXML(std::ostream&) const;
0348
0349     void SetRelErrorThreshold( double Threshold ) { fRelErrorThreshold = Threshold; }
0350     double GetRelErrorThreshold() const { return fRelErrorThreshold; }
0351
0352     void SetConstraintType( Constraint::Type ConstrType ) { fConstraintType = ConstrType; }
0353     Constraint::Type GetConstraintType() const { return fConstraintType; }
0354
0355   protected:
0356
0357     double fRelErrorThreshold;
0358     Constraint::Type fConstraintType;
0359
0360   };
0361
0362
0363 }
0364 }
0365
0366 #endif