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