]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - CORRFW/AliCFUnfolding.h
Fixing wrong usage of bits
[u/mrichter/AliRoot.git] / CORRFW / AliCFUnfolding.h
index bc611f6666e47b116302d550f23e93e5def9fd1e..0f2a3ad2ace9ede69ffed1b224db673ba2d46f2f 100644 (file)
@@ -14,6 +14,8 @@
 #include "THnSparse.h"
 
 class TF1;
+class TProfile;
+class TRandom3;
 
 class AliCFUnfolding : public TNamed {
 
@@ -25,9 +27,20 @@ class AliCFUnfolding : public TNamed {
   AliCFUnfolding& operator= (const AliCFUnfolding& c);
   ~AliCFUnfolding();
 
-  void SetMaxNumberOfIterations(Int_t n)  {fMaxNumIterations=n;}
-  void SetMaxChi2(Double_t val)           {fMaxChi2=val;}
-  void SetMaxChi2PerDOF(Double_t val);
+  void SetMaxNumberOfIterations(Int_t n = 10)  {fMaxNumIterations=n; fNRandomIterations=n; }
+
+  /* 
+     The following is for correct error estimation
+     thanks to Marta Verweij
+  */
+  void SetUseCorrelatedErrors  (Double_t maxConvergence = 1.e-06 , UInt_t randomSeed = 0)   {
+    fUseCorrelatedErrors = kTRUE             ;
+    fRandomSeed          = randomSeed        ;
+    SetMaxConvergencePerDOF(maxConvergence)  ;
+  }
+  void UnsetCorrelatedErrors() {fUseCorrelatedErrors=kFALSE;}
+
+
   void UseSmoothing(TF1* fcn=0x0, Option_t* opt="iremn") { // if fcn=0x0 then smooth using neighbouring bins 
     fUseSmoothing=kTRUE;                                   // this function must NOT be used if fNVariables > 3
     fSmoothFunction=fcn;                                   // the option "opt" is used if "fcn" is specified
@@ -36,16 +49,20 @@ class AliCFUnfolding : public TNamed {
                                                                                                 
   void Unfold();
 
-  THnSparse* GetResponse()        const {return fResponse;}
-  THnSparse* GetInverseResponse() const {return fInverseResponse;}
-  THnSparse* GetPrior()           const {return fPrior;}
-  THnSparse* GetOriginalPrior()   const {return fOriginalPrior;}
-  THnSparse* GetEfficiency()      const {return fEfficiency;}
-  THnSparse* GetUnfolded()        const {return fUnfolded;}
-  THnSparse* GetEstMeasured()     const {return fMeasuredEstimate;}
-  THnSparse* GetMeasured()        const {return fMeasured;}
-  THnSparse* GetConditional()     const {return fConditional;}
-  TF1*       GetSmoothFunction()  const {return fSmoothFunction;}
+  THnSparse* GetResponse()             const {return fResponse;}
+  THnSparse* GetInverseResponse()      const {return fInverseResponse;}
+  THnSparse* GetPrior()                const {return fPrior;}
+  THnSparse* GetOriginalPrior()        const {return fOriginalPrior;}
+  THnSparse* GetEfficiency()           const {return fEfficiency;}
+  THnSparse* GetUnfolded()             const {return fUnfolded;}
+  THnSparse* GetEstMeasured()          const {return fMeasuredEstimate;}
+  THnSparse* GetMeasured()             const {return fMeasured;}
+  THnSparse* GetConditional()          const {return fConditional;}
+  TF1*       GetSmoothFunction()       const {return fSmoothFunction;}
+  TProfile*  GetDeltaUnfoldedProfile() const {return fDeltaUnfoldedP;}
+  Int_t      GetDOF();                 // Returns number of degrees of freedom
+
+  static Short_t  SmoothUsingNeighbours(THnSparse*); // smoothes the unfolded spectrum using the neighbouring cells
 
  private :
   
@@ -58,14 +75,20 @@ class AliCFUnfolding : public TNamed {
                                       // then will be updated automatically at each iteration
   THnSparse     *fEfficiency;         // Efficiency map : dimensions must be N = number of variables
                                       // this map must be filled only with "true" values of the variables (should not include resolution effects)
-  THnSparse     *fMeasured;           // Measured spectrum to be unfolded : dimensions must be N = number of variables
+  THnSparse     *fMeasured;           // Measured spectrum to be unfolded : dimensions must be N = number of variables (modified)
+  THnSparse     *fMeasuredOrig;       // Measured spectrum to be unfolded : dimensions must be N = number of variables (not modified)
   Int_t          fMaxNumIterations;   // Maximum  number of iterations to be performed
   Int_t          fNVariables;         // Number of variables used in analysis spectra (pt, y, ...)
-  Double_t       fMaxChi2;            // Maximum Chi2 between unfolded and prior distributions. 
-  Bool_t         fUseSmoothing;       // Smooth the unfolded sectrum at each iteration
+/*   Double_t       fMaxChi2;            // Maximum Chi2 between unfolded and prior distributions.  */
+  Bool_t         fUseSmoothing;       // Smooth the unfolded sectrum at each iteration; default is kFALSE
   TF1           *fSmoothFunction;     // Function used to smooth the unfolded spectrum
   Option_t      *fSmoothOption;       // Option to use during the fit (with fSmoothFunction) ; default is "iremn"
 
+  /* correlated error calculation */
+  Double_t       fMaxConvergence;     // Convergence criterion in case of correlated error calculation
+  Bool_t         fUseCorrelatedErrors;// Calculate correlated errors for the final unfolded spectrum; default is kTRUE
+  Int_t          fNRandomIterations;  // Number of random distributed measured spectra
+
   // internal settings
   THnSparse     *fOriginalPrior;     // This is the original prior distribution : will not be modified
   THnSparse     *fInverseResponse;   // Inverse response matrix
@@ -76,7 +99,16 @@ class AliCFUnfolding : public TNamed {
   Int_t         *fCoordinates2N;     // Coordinates in 2N (measured,true) space
   Int_t         *fCoordinatesN_M;    // Coordinates in measured space
   Int_t         *fCoordinatesN_T;    // Coordinates in true space
-  
+
+
+  /* correlated error calculation */
+  THnSparse     *fRandomizedDist;    // Randomized distribution for each bin of the measured spectrum to calculate correlated errors
+  TRandom3      *fRandom3;           // Object to get random number following Poisson distribution
+  THnSparse     *fRandomUnfolded;
+  TProfile      *fDeltaUnfoldedP;    // Profile of the delta-unfolded distribution
+  Int_t          fNCalcCorrErrors;   // book keeping to prevend infinite loop
+  UInt_t         fRandomSeed;        // Random seed
+
 
   // functions
   void     Init();                  // initialisation of the internal settings
@@ -88,10 +120,17 @@ class AliCFUnfolding : public TNamed {
   void     CreateFlatPrior();       // creates a flat a priori distribution in case the one given in the constructor is null
   Double_t GetChi2();               // returns the chi2 between unfolded and prior spectra
   Short_t  Smooth();                // function calling smoothing methods
-  Short_t  SmoothUsingNeighbours(); // smoothes the unfolded spectrum using the neighbouring cells
   Short_t  SmoothUsingFunction();   // smoothes the unfolded spectrum using a fit function
 
-  ClassDef(AliCFUnfolding,0);
+  /* correlated error calculation */
+  Double_t GetConvergence();            // Returns convergence criterion
+  void     CalculateCorrelatedErrors(); // Calculates correlated errors for the final unfolded spectrum
+  void     InitDeltaUnfoldedProfile();  // Initializes the fDeltaUnfoldedP Profiles with spread option
+  void     CreateRandomizedDist();      // Create randomized dist from measured distribution
+  void     FillDeltaUnfoldedProfile();  // Fills the fDeltaUnfoldedP profile
+  void     SetMaxConvergencePerDOF (Double_t val);
+
+  ClassDef(AliCFUnfolding,1);
 };
 
 #endif