CORRFW/AliCFUnfolding.h

   1
   2 #ifndef ALICFUNFOLDING_H
   3 #define ALICFUNFOLDING_H
   4
   5 //--------------------------------------------------------------------//
   6 //                                                                    //
   7 // AliCFUnfolding Class                                               //
   8 // Class to handle general unfolding procedure using bayesian method  //
   9 //                                                                    //
  10 // Author : renaud.vernet@cern.ch                                     //
  11 //--------------------------------------------------------------------//
  12
  13 #include "TNamed.h"
  14 #include "THnSparse.h"
  15 #include "AliLog.h"
  16
  17 class TF1;
  18 class TRandom3;
  19
  20 class AliCFUnfolding : public TNamed {
  21
  22  public :
  23
  24   AliCFUnfolding();
  25   AliCFUnfolding(const Char_t* name, const Char_t* title, const Int_t nVar,
  26                  const THnSparse* response, const THnSparse* efficiency, const THnSparse* measured, const THnSparse* prior=0x0,
  27                  Double_t maxConvergencePerDOF = 1.e-06, UInt_t randomSeed = 0,
  28                  Int_t maxNumIterations = 10);
  29   ~AliCFUnfolding();
  30   void UnsetCorrelatedErrors()  {AliError("===================> DEPRECATED <=====================");}
  31   void SetUseCorrelatedErrors() {AliError("===================> DEPRECATED <=====================");}
  32   void SetMaxNumberOfIterations(Int_t n = 10)  {
  33     AliError("===================> DEPRECATED : should be set in constructor <=====================");
  34     AliError("===================> DEPRECATED : will be removed soon         <=====================");
  35     fMaxNumIterations = n;
  36   }
  37
  38   void SetNRandomIterations(Int_t n = 100) {fNRandomIterations = n;};
  39
  40   void UseSmoothing(TF1* fcn=0x0, Option_t* opt="iremn") { // if fcn=0x0 then smooth using neighbouring bins
  41     fUseSmoothing=kTRUE;                                   // this function must NOT be used if fNVariables > 3
  42     fSmoothFunction=fcn;                                   // the option "opt" is used if "fcn" is specified
  43     fSmoothOption=opt;
  44   }
  45
  46   void Unfold();
  47
  48   const THnSparse* GetResponse()             const {return fResponseOrig;}
  49   const THnSparse* GetEfficiency()           const {return fEfficiencyOrig;}
  50   const THnSparse* GetMeasured()             const {return fMeasuredOrig;}
  51   const THnSparse* GetOriginalPrior()        const {return fPriorOrig;}
  52         THnSparse* GetInverseResponse()      const {return fInverseResponse;}
  53         THnSparse* GetPrior()                const {return fPrior;}
  54         THnSparse* GetUnfolded()             const {return fUnfoldedFinal;}
  55         THnSparse* GetEstMeasured()          const {return fMeasuredEstimate;}
  56         THnSparse* GetConditional()          const {return fConditional;}
  57         TF1*       GetSmoothFunction()       const {return fSmoothFunction;}
  58         THnSparse* GetDeltaUnfoldedProfile() const {return fDeltaUnfoldedP;}
  59         Int_t      GetDOF();                 // Returns number of degrees of freedom
  60
  61   static Short_t  SmoothUsingNeighbours(THnSparse*); // smoothes the unfolded spectrum using the neighbouring cells
  62
  63  private :
  64   AliCFUnfolding(const AliCFUnfolding& c);
  65   AliCFUnfolding& operator= (const AliCFUnfolding& c);
  66
  67   //
  68   // user-related settings
  69   //
  70   const THnSparse     *fResponseOrig;     // Response matrix : dimensions must be 2N = 2 x (number of variables)
  71                                           // dimensions 0 ->  N-1 must be filled with reconstructed values
  72                                           // dimensions N -> 2N-1 must be filled with generated values
  73   const THnSparse     *fPriorOrig;        // This is the assumed generated distribution : dimensions must be N = number of variables
  74                                           // it will be used at the first step
  75                                           // then will be updated automatically at each iteration
  76   const THnSparse     *fEfficiencyOrig;   // Efficiency map : dimensions must be N = number of variables (modified)
  77                                           // this map must be filled only with "true" values of the variables (should not do "bin smearing")
  78   const THnSparse     *fMeasuredOrig;     // Measured spectrum to be unfolded : dimensions must be N = number of variables (modified)
  79
  80         Int_t          fMaxNumIterations; // Maximum  number of iterations to be performed
  81         Int_t          fNVariables;       // Number of variables used in analysis spectra (pt, y, ...)
  82         Bool_t         fUseSmoothing;     // Smooth the unfolded sectrum at each iteration; default is kFALSE
  83         TF1           *fSmoothFunction;   // Function used to smooth the unfolded spectrum
  84         Option_t      *fSmoothOption;     // Option to use during the fit (with fSmoothFunction) ; default is "iremn"
  85
  86   //
  87   // internal settings
  88   //
  89         Double_t       fMaxConvergence;    // Convergence criterion in case of correlated error calculation
  90         Int_t          fNRandomIterations; // Number of random distributed measured spectra
  91         THnSparse     *fResponse;          // Copy of the original response matrix    (modified)
  92         THnSparse     *fPrior;             // Copy of the original prior spectrum     (modified)
  93         THnSparse     *fEfficiency;        // Copy of original efficiency             (modified)
  94         THnSparse     *fMeasured;          // Copy of the original measureed spectrum (modified)
  95         THnSparse     *fInverseResponse;   // Inverse response matrix
  96         THnSparse     *fMeasuredEstimate;  // Estimation of the measured (M) spectrum given the a priori (T) distribution
  97         THnSparse     *fConditional;       // Matrix holding the conditional probabilities P(M|T)
  98         THnSparse     *fUnfolded;          // Unfolded spectrum (modified before and during error calculation)
  99         THnSparse     *fUnfoldedFinal;     // Final unfolded spectrum
 100         Int_t         *fCoordinates2N;     // Coordinates in 2N (measured,true) space
 101         Int_t         *fCoordinatesN_M;    // Coordinates in measured space
 102         Int_t         *fCoordinatesN_T;    // Coordinates in true space
 103
 104
 105   /* correlated error calculation */
 106   THnSparse     *fRandomResponse;    // Randomized distribution for each bin of the response matrix     to calculate correlated errors
 107   THnSparse     *fRandomEfficiency;  // Randomized distribution for each bin of the efficiency spectrum to calculate correlated errors
 108   THnSparse     *fRandomMeasured;    // Randomized distribution for each bin of the measured   spectrum to calculate correlated errors
 109   TRandom3      *fRandom3;           // Object to get random number following Poisson distribution
 110   THnSparse     *fDeltaUnfoldedP;    // Profile of the delta-unfolded distribution
 111   THnSparse     *fDeltaUnfoldedN;    // Entries of the delta-unfolded distribution (count for each bin)
 112   Short_t        fNCalcCorrErrors;   // Book-keeping to prevend infinite loop
 113   UInt_t         fRandomSeed;        // Random seed
 114
 115
 116   // functions
 117   void     Init();                  // initialisation of the internal settings
 118   void     GetCoordinates();        // gets a cell coordinates in Measured and True space
 119   void     CreateConditional();     // creates the conditional matrix from the response matrix
 120   void     CreateEstMeasured();     // creates the measured spectrum estimation from the conditional matrix and the prior distribution
 121   void     CreateInvResponse();     // creates the inverse response function (Bayes Theorem) from the conditional matrix and the prior distribution
 122   void     CreateUnfolded();        // creates the unfolded spectrum from the inverse response matrix and the measured distribution
 123   void     CreateFlatPrior();       // creates a flat a priori distribution in case the one given in the constructor is null
 124   Double_t GetChi2();               // returns the chi2 between unfolded and prior spectra
 125   Short_t  Smooth();                // function calling smoothing methods
 126   Short_t  SmoothUsingFunction();   // smoothes the unfolded spectrum using a fit function
 127
 128   /* correlated error calculation */
 129   Double_t GetConvergence();            // Returns convergence criterion
 130   void     CalculateCorrelatedErrors(); // Calculates correlated errors for the final unfolded spectrum
 131   void     CreateRandomizedDist();      // Create randomized dist from measured distribution
 132   void     FillDeltaUnfoldedProfile();  // Fills the fDeltaUnfoldedP profile
 133   void     SetMaxConvergencePerDOF (Double_t val);
 134
 135   ClassDef(AliCFUnfolding,1);
 136 };
 137
 138 #endif