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
0208     ShapeSys() :
0209       HistogramUncertaintyBase(),
0210       fConstraintType(Constraint::Gaussian) {}
0211     ShapeSys(const ShapeSys& other) :
0212       HistogramUncertaintyBase(other),
0213       fConstraintType(other.fConstraintType) {}
0214     ShapeSys& operator=(const ShapeSys& oth) {
0215        if (this == &oth) return *this;
0216        HistogramUncertaintyBase::operator=(oth);
0217        fConstraintType = oth.fConstraintType;
0218        return *this;
0219     }
0220     ShapeSys& operator=(ShapeSys&&) = default;
0221
0222     void SetInputFile( const std::string& InputFile ) { fInputFileHigh = InputFile; }
0223     std::string GetInputFile() const { return fInputFileHigh; }
0224
0225     void SetHistoName( const std::string& HistoName ) { fHistoNameHigh = HistoName; }
0226     std::string GetHistoName() const { return fHistoNameHigh; }
0227
0228     void SetHistoPath( const std::string& HistoPath ) { fHistoPathHigh = HistoPath; }
0229     std::string GetHistoPath() const { return fHistoPathHigh; }
0230
0231     void Print(std::ostream& = std::cout) const override;
0232     void PrintXML(std::ostream&) const override;
0233     void writeToFile( const std::string& FileName, const std::string& DirName ) override;
0234
0235     const TH1* GetErrorHist() const {
0236       return fhHigh.get();
0237     }
0238     void SetErrorHist(TH1* hError) {
0239       fhHigh.reset(hError);
0240     }
0241
0242     void SetConstraintType( Constraint::Type ConstrType ) { fConstraintType = ConstrType; }
0243     Constraint::Type GetConstraintType() const { return fConstraintType; }
0244
0245   protected:
0246     Constraint::Type fConstraintType;
0247   };
0248
0249 /** \class ShapeFactor
0250  * \ingroup HistFactory
0251  * *Un*constrained bin-by-bin variation of affected histogram.
0252  */
0253   class ShapeFactor : public HistogramUncertaintyBase {
0254
0255   public:
0256
0257     ShapeFactor() :
0258       HistogramUncertaintyBase(),
0259       fConstant{false},
0260       fHasInitialShape{false} {}
0261
0262     void Print(std::ostream& = std::cout) const override;
0263     void PrintXML(std::ostream&) const override;
0264     void writeToFile( const std::string& FileName, const std::string& DirName) override;
0265
0266     void SetInitialShape(TH1* shape) {
0267       fhHigh.reset(shape);
0268     }
0269     const TH1* GetInitialShape() const { return fhHigh.get(); }
0270
0271     void SetConstant(bool constant) { fConstant = constant; }
0272     bool IsConstant() const { return fConstant; }
0273
0274     bool HasInitialShape() const { return fHasInitialShape; }
0275
0276     void SetInputFile( const std::string& InputFile ) {
0277       fInputFileHigh = InputFile;
0278       fHasInitialShape=true;
0279     }
0280     const std::string& GetInputFile() const { return fInputFileHigh; }
0281
0282     void SetHistoName( const std::string& HistoName ) {
0283       fHistoNameHigh = HistoName;
0284       fHasInitialShape=true;
0285     }
0286     const std::string& GetHistoName() const { return fHistoNameHigh; }
0287
0288     void SetHistoPath( const std::string& HistoPath ) {
0289       fHistoPathHigh = HistoPath;
0290       fHasInitialShape=true;
0291     }
0292     const std::string& GetHistoPath() const { return fHistoPathHigh; }
0293
0294   protected:
0295
0296     bool fConstant;
0297
0298     // A histogram representing
0299     // the initial shape
0300     bool fHasInitialShape;
0301   };
0302
0303 /** \class StatError
0304  * \ingroup HistFactory
0305  * Statistical error of Monte Carlo predictions.
0306  */
0307   class StatError : public HistogramUncertaintyBase {
0308
0309   public:
0310
0311     StatError() :
0312       HistogramUncertaintyBase(),
0313       fActivate(false), fUseHisto(false) {}
0314
0315     void Print(std::ostream& = std::cout) const override;
0316     void PrintXML(std::ostream&) const override;
0317     void writeToFile( const std::string& FileName, const std::string& DirName ) override;
0318
0319     void Activate( bool IsActive=true ) { fActivate = IsActive; }
0320     bool GetActivate() const { return fActivate; }
0321
0322     void SetUseHisto( bool UseHisto=true ) { fUseHisto = UseHisto; }
0323     bool GetUseHisto() const { return fUseHisto; }
0324
0325     void SetInputFile( const std::string& InputFile ) { fInputFileHigh = InputFile; }
0326     const std::string& GetInputFile() const { return fInputFileHigh; }
0327
0328     void SetHistoName( const std::string& HistoName ) { fHistoNameHigh = HistoName; }
0329     const std::string& GetHistoName() const { return fHistoNameHigh; }
0330
0331     void SetHistoPath( const std::string& HistoPath ) { fHistoPathHigh = HistoPath; }
0332     const std::string& GetHistoPath() const { return fHistoPathHigh; }
0333
0334
0335     const TH1* GetErrorHist() const {
0336       return fhHigh.get();
0337     }
0338     void SetErrorHist(TH1* Error) {
0339       fhHigh.reset(Error);
0340     }
0341
0342   protected:
0343
0344     bool fActivate;
0345     bool fUseHisto; // Use an external histogram for the errors
0346   };
0347
0348 /** \class StatErrorConfig
0349  * \ingroup HistFactory
0350  * Configuration to automatically assign nuisance parameters for the statistical
0351  * error of the Monte Carlo simulations.
0352  * The default is to assign a Poisson uncertainty to a bin when its statistical uncertainty
0353  * is larger than 5% of the bin content.
0354  */
0355   class StatErrorConfig {
0356
0357   public:
0358
0359     StatErrorConfig() : fRelErrorThreshold( .05 ), fConstraintType( Constraint::Poisson ) {;}
0360     void Print(std::ostream& = std::cout) const;
0361     void PrintXML(std::ostream&) const;
0362
0363     void SetRelErrorThreshold( double Threshold ) { fRelErrorThreshold = Threshold; }
0364     double GetRelErrorThreshold() const { return fRelErrorThreshold; }
0365
0366     void SetConstraintType( Constraint::Type ConstrType ) { fConstraintType = ConstrType; }
0367     Constraint::Type GetConstraintType() const { return fConstraintType; }
0368
0369   protected:
0370
0371     double fRelErrorThreshold;
0372     Constraint::Type fConstraintType;
0373
0374   };
0375
0376
0377 }
0378 }
0379
0380 #endif