/cvmfs/atlas.cern.ch/repo/sw/ASG/AnalysisBase/2.4.31/ElectronPhotonSelectorTools/ElectronPhotonSelectorTools/TElectronLikelihoodTool.h

// Dear emacs, this is -*-c++-*-
//----------------------------------------------------------------------------------------

//-------------------------------------------------------------------------------------------------------


#ifndef TELECTRONLIKELIHOODTOOL_H
#define TELECTRONLIKELIHOODTOOL_H

#include <fstream>
#include <iostream>
// Include the return objects and the base class
#include "PATCore/TCalculatorToolBase.h"
#include "PATCore/TSelectorToolBase.h"
#include "AsgTools/AsgMessaging.h"

#include <string>                       // for string
#include <vector>                       // for vector
class TFile;
class TH1F;
namespace Root { class TAccept; }
namespace Root { class TResult; }
namespace{
  const unsigned int  IP_BINS=1;
}
namespace LikeEnum {

  struct LHAcceptVars_t{
    double likelihood;
    double eta;
    double eT;
    int nSiHitsPlusDeadSensors;
    int nPixHitsPlusDeadSensors;
    bool passBLayerRequirement;
    int convBit;
    double d0;
    double deltaEta;
    double deltaphires;
    double wstot;
    double EoverP;
    double ip;         
  };

  struct LHCalcVars_t{
    double eta;
    double eT;
    double f3; 
    double rHad; 
    double rHad1;
    double Reta;
    double w2; 
    double f1; 
    double eratio;
    double deltaEta; 
    double d0; 
    double d0sigma;
    double rphi; 
    double deltaPoverP;
    double deltaphires;
    double TRT_PID;
    double ip;        
  };
}

namespace Root {
  class TElectronLikelihoodTool : public TCalculatorToolBase, public TSelectorToolBase, public asg::AsgMessaging
  {

  public: 

    TElectronLikelihoodTool(const char* name = "TElectronLikelihoodTool");
    
    ~TElectronLikelihoodTool();

  private:
    class SafeTH1{
    public :
      SafeTH1(TH1F* hist);
      ~SafeTH1();
      int GetNbinsX();
      int FindBin(double);
      double GetBinContent(int);
      double GetBinLowEdge(int);
      double Integral();

    private:
      std::vector<float> m_binContent;
      double m_firstBinLowEdge;
      double m_lastBinLowEdge;
      double m_binWidth;
      double m_integral;
    };

    // Main methods
  public:
    int initialize();

    int finalize() {return 1;}

    const Root::TAccept& accept(LikeEnum::LHAcceptVars_t& vars_struct) const;
    const Root::TAccept& accept( double likelihood,
                                 double eta, double eT,
                                 int nSiHitsPlusDeadSensors, int nPixHitsPlusDeadSensors,
                                 bool passBLayerRequirement,
                                 int convBit, double d0, double deltaEta, double deltaphires, 
                                 double wstot, double EoverP, double ip ) const;
    const Root::TResult& calculate(LikeEnum::LHCalcVars_t& vars_struct) const ;
    const Root::TResult& calculate( double eta, double eT,double f3, double rHad, double rHad1,
                                    double Reta, double w2, double f1, double eratio,
                                    double deltaEta, double d0, double d0sigma, double rphi,
                                    double deltaPoverP ,double deltaphires, double TRT_PID,
                                    double ip) const;
    

    inline void setPDFFileName ( const std::string& val ) { PdfFileName = val; }

    inline void setVariableNames ( const std::string& val ) { 
      VariableNames = val; 
      m_variableBitMask = GetLikelihoodBitmask(val);
    }

    int LoadVarHistograms(std::string vstr, unsigned int varIndex);

    inline void setBinning ( const std::string& val ) { m_ipBinning = val; }

    inline void setResultPrefix ( const std::string& val ) { m_resultPrefix = val; }
   
    inline void setResultName ( const std::string& val ) { m_resultName = val; }

    unsigned int getBitmask(void) const { return m_variableBitMask;} 
    inline void setBitmask(unsigned int val) { m_variableBitMask = val; };

    // Private methods
  private:
   
    // For every input "varVector", make sure elements of vector are
    // in the same order as prescribed in fVariables

    double EvaluateLikelihood(std::vector<double> varVector,double et,double eta,double ip=0) const;

    double EvaluateLikelihood(std::vector<float>  varVector,double et,double eta,double ip=0) const;


    // To concoct a bitmask on your own, use the 
    // variable names prescribed in fVariables.
    
    unsigned int GetLikelihoodBitmask(std::string vars) const;
    
    double InterpolateCuts(const std::vector<double>& cuts,const std::vector<double>& cuts_4gev,double et,double eta) const;
    double InterpolatePdfs(unsigned int s_or_b,unsigned int ipbin,double et,double eta,int bin,unsigned int var) const;
    
  public:
    std::vector<int> CutBL;
    std::vector<int> CutPi;
    std::vector<int> CutSi;
    std::vector<double> CutA0;
    std::vector<double> CutDeltaEta;
    // /** @brief do cut on delta phi bit*/
    std::vector<double> CutDeltaPhiRes;
    bool doCutConversion;
    bool doRemoveF3AtHighEt;
    bool doRemoveTRTPIDAtHighEt;
    bool doSmoothBinInterpolation;
    bool useHighETLHBinning;
    bool useOneExtraHighETLHBin;
    double HighETBinThreshold;
    // /** @brief do cut on wstot above HighETBinThreshold bit*/
    std::vector<double> CutWstotAtHighET;
    // /** @brief do cut on EoverP above HighETBinThreshold bit*/
    std::vector<double> CutEoverPAtHighET;
    bool doPileupTransform;
    bool doCentralityTransform;
    std::vector<double> CutLikelihood;
    std::vector<double> CutLikelihoodPileupCorrection;
    std::vector<double> CutLikelihood4GeV;
    std::vector<double> CutLikelihoodPileupCorrection4GeV;
    std::vector<double> DiscHardCutForPileupTransform;
    std::vector<double> DiscHardCutSlopeForPileupTransform;
    std::vector<double> DiscHardCutQuadForPileupTransform;
    std::vector<double> DiscLooseForPileupTransform;
    std::vector<double> DiscHardCutForPileupTransform4GeV;
    std::vector<double> DiscHardCutSlopeForPileupTransform4GeV;
    std::vector<double> DiscHardCutQuadForPileupTransform4GeV;
    std::vector<double> DiscLooseForPileupTransform4GeV;
    double DiscMaxForPileupTransform;
    double PileupMaxForPileupTransform;
    std::string VariableNames;
    std::string PdfFileName;


    // Private methods
  private:
    double TransformLikelihoodOutput(double ps,double pb, double ip, double et, double eta) const;

    unsigned int getLikelihoodEtaBin(double eta) const ;

    unsigned int getLikelihoodEtHistBin(double eT)const ;
    
    unsigned int getLikelihoodEtDiscBin(double eT) const;


    // Private member variables
  private:

    unsigned int        m_variableBitMask;

    std::string         m_ipBinning;

    TFile*              m_pdfFile;

    std::string m_resultPrefix;

    std::string m_resultName;

    int m_cutPosition_kinematic;

    int m_cutPosition_NSilicon;

    int m_cutPosition_NPixel;

    int m_cutPosition_NBlayer;

    int m_cutPosition_conversion;

    int m_cutPosition_LH;

    int m_cutPositionTrackA0;

    int m_cutPositionTrackMatchEta;

    // /// The position of the deltaphi cut bit in the TAccept return object
    int m_cutPositionTrackMatchPhiRes;

    // /// The position of the high ET wstot cut bit in the TAccept return object
    int m_cutPositionWstotAtHighET;

    // /// The position of the high ET EoverP cut bit in the TAccept return object
    int m_cutPositionEoverPAtHighET;

    int m_resultPosition_LH;

    static const double fIpBounds[IP_BINS+1];
    static const unsigned int  fnEtBinsHist     = 8;  // number of hists stored for LH with many high ET bins (useHighETLHBinning), including 4GeV bin
    static const unsigned int  fnDiscEtBins     = 33; // number of discs stored for LH with many high ET bins (useHighETLHBinning), excluding 4GeV bin
    static const unsigned int  fnEtBinsHistOrig = 7;  // number of hists stored for original LH, including 4GeV bin (for backwards compatibility)
    static const unsigned int  fnDiscEtBinsOrig = 9;  // number of discs stored for original LH, excluding 4GeV bin (for backwards compatibility)
    static const unsigned int  fnDiscEtBinsOneExtra = 10; // number of discs stored for original LH plus one for HighETBinThreshold (useOneExtraHighETLHBin), excluding 4GeV bin
    static const unsigned int  fnEtaBins        = 10;
    static const unsigned int  fnVariables      = 13;
    TElectronLikelihoodTool::SafeTH1*      fPDFbins     [2][IP_BINS][fnEtBinsHist][fnEtaBins][fnVariables]; // [sig(0)/bkg(1)][ip][et][eta][variable]
    static const std::string  fVariables                [fnVariables];

    unsigned int getIpBin(double ip) const;
    void getBinName(char* buffer, int etbin,int etabin, int ipbin, std::string iptype) const;
  };

} // End: namespace Root

//----------------------------------------------------------------------------------------
#endif

---