PWG4/JetTasks/AliJetSpectrumUnfolding.h

   1 #ifndef ALIJETSPECTRUMUNFOLDING_H
   2 #define ALIJETSPECTRUMUNFOLDING_H
   3
   4 #include "TNamed.h"
   5 //
   6 // class that contains the correction matrix and the functions for
   7 // correction the jet spectrum
   8 // implements 2-dim bayesian method
   9 //
  10
  11 class TH1;
  12 class TH2;
  13 class TH3;
  14 class TH1F;
  15 class TH2F;
  16 class TH3F;
  17 class TF1;
  18 class TF2;
  19 class TCollection;
  20
  21 #include <TMatrixD.h>
  22 #include <TVectorD.h>
  23 #include <THnSparse.h>
  24
  25 class AliJetSpectrumUnfolding : public TNamed {
  26   public:
  27     AliJetSpectrumUnfolding();
  28     AliJetSpectrumUnfolding(const Char_t* name, const Char_t* title);
  29     virtual ~AliJetSpectrumUnfolding();
  30
  31     virtual Long64_t Merge(TCollection* list);
  32
  33     Bool_t LoadHistograms(const Char_t* dir = 0);
  34     void SaveHistograms();
  35     void DrawHistograms();
  36     void DrawComparison(const char* name, TH2* genSpec);
  37
  38     void SetBayesianParameters(Float_t smoothing, Int_t nIterations);
  39     void ApplyBayesianMethod(Float_t regPar = 1, Int_t nIterations = 100, TH2* initialConditions = 0, Bool_t determineError = kTRUE);
  40
  41     TH2F* GetRecSpectrum()      { return fRecSpectrum; }
  42     TH2F* GetGenSpectrum()      { return fGenSpectrum; }
  43     TH2F* GetUnfSpectrum()      { return fUnfSpectrum; }
  44     THnSparseF* GetCorrelation(){ return fCorrelation; }
  45
  46     void SetRecSpectrum(TH2F* hist)      { fRecSpectrum  = hist; }
  47     void SetGenSpectrum(TH2F* hist)      { fGenSpectrum  = hist; }
  48     void SetUnfSpectrum(TH2F* hist)      { fUnfSpectrum  = hist; }
  49     void SetCorrelation(THnSparseF* hist){ fCorrelation  = hist; }
  50
  51     void SetGenRecFromFunc(TF2* inputGen);
  52     TH2F* CalculateRecSpectrum(TH2* inputGen);
  53
  54     static void NormalizeToBinWidth(TH2* hist);
  55
  56     void GetComparisonResults(Float_t* gen = 0, Int_t* genLimit = 0, Float_t* residuals = 0, Float_t* ratioAverage = 0) const;
  57
  58   protected:
  59     void SetupCurrentHists(Bool_t createBigBin);
  60
  61     static Double_t BayesCov(THnSparseF* M, THnSparseF* correlation, Int_t* binTM, Int_t* binTM1);
  62     static Double_t BayesUncertaintyTerms(THnSparseF *M, THnSparseF *C, Int_t* binTM, Int_t* binTM1, Double_t Nt);
  63     static Int_t UnfoldWithBayesian(THnSparseF* correlation, TH2* measured, TH2* initialConditions, TH2* aResult, Float_t regPar, Int_t nIterations, Bool_t calculateErrors = kFALSE);
  64
  65     static Float_t fgBayesianSmoothing;             //! smoothing parameter (0 = no smoothing)
  66     static Int_t   fgBayesianIterations;            //! number of iterations in Bayesian method
  67
  68     TH2F* fCurrentRec;
  69     THnSparseF* fCurrentCorrelation;
  70
  71     TH2F* fRecSpectrum;
  72     TH2F* fGenSpectrum;
  73     TH2F* fUnfSpectrum;
  74     THnSparseF* fCorrelation;
  75
  76     Float_t fLastChi2MC;        //! last Chi2 between MC and unfolded ESD (calculated in DrawComparison)
  77     Int_t   fLastChi2MCLimit;   //! bin where the last chi2 breached a certain threshold
  78     Float_t fLastChi2Residuals; //! last Chi2 of the ESD and the folded unfolded ESD (calculated in DrawComparison)
  79     Float_t fRatioAverage;      //! last average of |ratio-1| where ratio = unfolded / mc (bin 2..150)
  80
  81  private:
  82     AliJetSpectrumUnfolding(const AliJetSpectrumUnfolding&);
  83     AliJetSpectrumUnfolding& operator=(const AliJetSpectrumUnfolding&);
  84
  85   ClassDef(AliJetSpectrumUnfolding, 2);
  86 };
  87
  88 #endif
  89