e6339065 |
1 | /* $Id$ */ |
2 | |
3 | #ifndef ALIUNFOLDING_H |
4 | #define ALIUNFOLDING_H |
5 | |
6 | // |
7 | // class that implements several unfolding methods |
8 | // E.g. chi2 minimization and bayesian unfolding |
9 | // |
10 | |
11 | // TMatrixD, TVectorD defined here, because it does not seem possible to predeclare these (or i do not know how) |
12 | // --> |
13 | // $ROOTSYS/include/TVectorDfwd.h:21: conflicting types for `typedef struct TVectorT<Double_t> TVectorD' |
14 | // PWG0/AliUnfolding.h:21: previous declaration as `struct TVectorD' |
15 | |
16 | #include "TObject.h" |
17 | #include <TMatrixD.h> |
18 | #include <TVectorD.h> |
19 | |
20 | class TH1; |
21 | class TH2; |
22 | class TH1F; |
23 | class TH2F; |
24 | class TH3F; |
25 | class TF1; |
26 | class TCollection; |
27 | |
28 | class AliUnfolding : public TObject |
29 | { |
30 | public: |
31 | enum RegularizationType { kNone = 0, kPol0, kPol1, kLog, kEntropy, kCurvature }; |
32 | enum MethodType { kChi2Minimization = 0, kBayesian = 1 }; |
33 | |
34 | AliUnfolding(); |
35 | virtual ~AliUnfolding(); |
36 | |
37 | void SetInput(TH2* correlationMatrix, TH1* efficiency, TH1* measured) { fCurrentCorrelation = correlationMatrix; fCurrentEfficiency = efficiency; fCurrentESD = measured; } |
38 | void SetInitialConditions(TH1* initialConditions) { fInitialConditions = initialConditions; } |
39 | const TH1* GetResult() const { return fResult; } |
40 | |
41 | static void SetParameters(Int_t measuredBins, Int_t unfoldedBins, Bool_t bigbin) { fMaxInput = measuredBins; fMaxParams = unfoldedBins; fgCreateBigBin = bigbin; } |
42 | static void SetChi2MinimizationParameters(RegularizationType type, Float_t weight) { fgRegularizationType = type; fgRegularizationWeight = weight; } |
43 | static void SetRegularizationRange(Int_t start, Int_t end) { fgRegularizationRangeStart = start; fgRegularizationRangeEnd = end; } |
44 | static void SetBayesianParameters(Float_t smoothing, Int_t nIterations) { fgBayesianSmoothing = smoothing; fgBayesianIterations = nIterations; } |
45 | |
46 | Int_t ApplyMinuitFit(Bool_t check = kFALSE); |
47 | Int_t ApplyBayesianMethod(Bool_t determineError = kTRUE); |
48 | Int_t ApplyNBDFit(); |
49 | Int_t ApplyLaszloMethod(); |
50 | |
51 | TH1* StatisticalUncertainty(MethodType methodType, Bool_t randomizeMeasured, Bool_t randomizeResponse, TH1* compareTo = 0); |
52 | |
53 | protected: |
54 | static Double_t RegularizationPol0(TVectorD& params); |
55 | static Double_t RegularizationPol1(TVectorD& params); |
56 | static Double_t RegularizationTotalCurvature(TVectorD& params); |
57 | static Double_t RegularizationEntropy(TVectorD& params); |
58 | static Double_t RegularizationLog(TVectorD& params); |
59 | |
60 | static void MinuitFitFunction(Int_t&, Double_t*, Double_t& chi2, Double_t *params, Int_t); |
61 | static void MinuitNBD(Int_t& unused1, Double_t* unused2, Double_t& chi2, Double_t *params, Int_t unused3); |
62 | |
63 | void SetupCurrentHists(); |
64 | |
65 | Int_t UnfoldWithBayesian(* aEfficiency, TH1* measured, TH1* initialConditions, TH1* aResult, Float_t regPar, Int_t nIterations); |
66 | Int_t UnfoldWithMinuit(TH1* correlation, TH1* aEfficiency, TH1* measured, TH1* initialConditions, TH1* result, Bool_t check); |
67 | |
68 | Float_t BayesCovarianceDerivate(Float_t matrixM[251][251], TH2* hResponse, Int_t k, Int_t i, Int_t r, Int_t u); |
69 | |
70 | TH1* fCurrentESD; //! measured spectrum |
71 | TH2* fCurrentCorrelation; //! correlation matrix |
72 | TH1* fCurrentEfficiency; //! efficiency |
73 | TH1* fInitialConditions; //! initial conditions |
74 | TH1* fResult; //! unfolding result |
75 | |
76 | // static variable to be accessed by MINUIT |
77 | static TMatrixD* fgCorrelationMatrix; //! contains fCurrentCorrelation in matrix form |
78 | static TMatrixD* fgCorrelationCovarianceMatrix; //! contains the errors of fCurrentESD |
79 | static TVectorD* fgCurrentESDVector; //! contains fCurrentESD |
80 | static TVectorD* fgEntropyAPriori; //! a-priori distribution for entropy regularization |
81 | |
82 | static TF1* fgNBD; //! negative binomial distribution |
83 | |
84 | static Int_t fgMaxParams; //! bins in unfolded histogram = number of fit params |
85 | static Int_t fgMaxInput; //! bins in measured histogram |
86 | |
87 | // configuration params follow |
88 | static RegularizationType fgRegularizationType; //! regularization that is used during Chi2 method |
89 | static Float_t fgRegularizationWeight; //! factor for regularization term |
90 | static Int_t fgRegularizationRangeStart; //! first bin where regularization is applied |
91 | static Int_t fgRegularizationRangeEnd; //! last bin + 1 where regularization is applied |
92 | static Bool_t fgCreateBigBin; //! to fix fluctuations at high multiplicities, all entries above a certain limit are summarized in one bin |
93 | |
94 | static Float_t fgBayesianSmoothing; //! smoothing parameter (0 = no smoothing) |
95 | static Int_t fgBayesianIterations; //! number of iterations in Bayesian method |
96 | // end of configuration |
97 | |
98 | private: |
99 | AliUnfolding(const AliUnfolding&); |
100 | AliUnfolding& operator=(const AliUnfolding&); |
101 | |
102 | ClassDef(AliUnfolding, 0); |
103 | }; |
104 | |
105 | #endif |
106 | |