]>
Commit | Line | Data |
---|---|---|
e6339065 | 1 | /* $Id$ */ |
2 | ||
3 | #ifndef ALIUNFOLDING_H | |
4 | #define ALIUNFOLDING_H | |
5 | ||
6 | // | |
7 | // class that implements several unfolding methods | |
19442b86 | 8 | // I.e. chi2 minimization and bayesian unfolding |
9 | // The whole class is static and not thread-safe (due to the fact that MINUIT unfolding is not thread-safe) | |
e6339065 | 10 | // |
11 | ||
12 | // TMatrixD, TVectorD defined here, because it does not seem possible to predeclare these (or i do not know how) | |
13 | // --> | |
14 | // $ROOTSYS/include/TVectorDfwd.h:21: conflicting types for `typedef struct TVectorT<Double_t> TVectorD' | |
15 | // PWG0/AliUnfolding.h:21: previous declaration as `struct TVectorD' | |
16 | ||
17 | #include "TObject.h" | |
18 | #include <TMatrixD.h> | |
19 | #include <TVectorD.h> | |
20 | ||
21 | class TH1; | |
22 | class TH2; | |
e6339065 | 23 | class TF1; |
9e065ad2 | 24 | class TCanvas; |
25 | class TVirtualPad; | |
26 | class TAxis; | |
e6339065 | 27 | |
28 | class AliUnfolding : public TObject | |
29 | { | |
30 | public: | |
9e065ad2 | 31 | enum RegularizationType { kNone = 0, kPol0, kPol1, kLog, kEntropy, kCurvature, kRatio, kPowerLaw, kLogLog }; |
19442b86 | 32 | enum MethodType { kInvalid = -1, kChi2Minimization = 0, kBayesian = 1, kFunction = 2}; |
e6339065 | 33 | |
19442b86 | 34 | virtual ~AliUnfolding() {}; |
e6339065 | 35 | |
19442b86 | 36 | static void SetUnfoldingMethod(MethodType methodType); |
37 | static void SetCreateOverflowBin(Float_t overflowBinLimit); | |
38 | static void SetSkipBinsBegin(Int_t nBins); | |
39 | static void SetNbins(Int_t nMeasured, Int_t nUnfolded); | |
40 | static void SetChi2Regularization(RegularizationType type, Float_t weight); | |
41 | static void SetMinuitStepSize(Float_t stepSize) { fgMinuitStepSize = stepSize; } | |
9e065ad2 | 42 | static void SetMinuitPrecision(Float_t pres) {fgMinuitPrecision = pres;} |
d22beb7f | 43 | static void SetMinuitMaxIterations(Int_t iter) {fgMinuitMaxIterations = iter;} |
9f070ef5 | 44 | static void SetMinuitStrategy(Double_t strat) {fgMinuitStrategy = strat;} |
95e970ca | 45 | static void SetMinimumInitialValue(Bool_t flag, Float_t value = -1) { fgMinimumInitialValue = flag; fgMinimumInitialValueFix = value; } |
19442b86 | 46 | static void SetNormalizeInput(Bool_t flag) { fgNormalizeInput = flag; } |
95e970ca | 47 | static void SetNotFoundEvents(Float_t notFoundEvents) { fgNotFoundEvents = notFoundEvents; } |
48 | static void SetSkip0BinInChi2(Bool_t flag) { fgSkipBin0InChi2 = flag; } | |
19442b86 | 49 | static void SetBayesianParameters(Float_t smoothing, Int_t nIterations); |
50 | static void SetFunction(TF1* function); | |
51 | static void SetDebug(Bool_t flag) { fgDebug = flag; } | |
9e065ad2 | 52 | static void SetPowern(Int_t n) {fgPowern = -1*n;} |
e6339065 | 53 | |
19442b86 | 54 | static Int_t Unfold(TH2* correlation, TH1* efficiency, TH1* measured, TH1* initialConditions, TH1* result, Bool_t check = kFALSE); |
95e970ca | 55 | static Int_t UnfoldGetBias(TH2* correlation, TH1* efficiency, TH1* measured, TH1* initialConditions, TH1* result); |
e6339065 | 56 | |
95e970ca | 57 | static TH1* GetPenaltyPlot(Double_t* params); |
58 | static TH1* GetPenaltyPlot(TH1* corrected); | |
9e065ad2 | 59 | |
60 | static TH1* GetResidualsPlot(Double_t* params); | |
61 | static TH1* GetResidualsPlot(TH1* corrected); | |
62 | ||
63 | static Double_t GetChi2FromFit() {return fChi2FromFit;} | |
64 | static Double_t GetPenaltyVal() {return fPenaltyVal;} | |
65 | static Double_t GetAvgResidual() {return fAvgResidual;} | |
66 | ||
67 | static void DrawResults(TH2* correlation, TH1* efficiency, TH1* measured, TH1* initialConditions, TCanvas *canvas = 0, Int_t reuseHists = kFALSE,TH1 *unfolded=0); | |
68 | static void InteractiveUnfold(TH2* correlation, TH1* efficiency, TH1* measured, TH1* initialConditions); | |
69 | static void RedrawInteractive(); | |
70 | ||
19442b86 | 71 | protected: |
72 | AliUnfolding() {}; | |
e6339065 | 73 | |
19442b86 | 74 | static Int_t UnfoldWithMinuit(TH2* correlation, TH1* efficiency, TH1* measured, TH1* initialConditions, TH1* result, Bool_t check); |
75 | static Int_t UnfoldWithBayesian(TH2* correlation, TH1* aEfficiency, TH1* measured, TH1* initialConditions, TH1* aResult); | |
76 | static Int_t UnfoldWithFunction(TH2* correlation, TH1* efficiency, TH1* measured, TH1* /* initialConditions */, TH1* aResult); | |
77 | ||
19442b86 | 78 | static void CreateOverflowBin(TH2* correlation, TH1* measured); |
95e970ca | 79 | static void SetStaticVariables(TH2* correlation, TH1* measured, TH1* efficiency); |
e6339065 | 80 | |
9e065ad2 | 81 | static void MakePads(); |
82 | static void DrawGuess(Double_t *params, TVirtualPad *pfolded=0, TVirtualPad *pres=0, TVirtualPad *ppen=0, Int_t reuseHists = kFALSE, TH1* unfolded=0); | |
83 | ||
e6339065 | 84 | static Double_t RegularizationPol0(TVectorD& params); |
85 | static Double_t RegularizationPol1(TVectorD& params); | |
86 | static Double_t RegularizationTotalCurvature(TVectorD& params); | |
87 | static Double_t RegularizationEntropy(TVectorD& params); | |
88 | static Double_t RegularizationLog(TVectorD& params); | |
95e970ca | 89 | static Double_t RegularizationRatio(TVectorD& params); |
9e065ad2 | 90 | static Double_t RegularizationPowerLaw(TVectorD& params); |
91 | static Double_t RegularizationLogLog(TVectorD& params); | |
e6339065 | 92 | |
19442b86 | 93 | static void Chi2Function(Int_t&, Double_t*, Double_t& chi2, Double_t *params, Int_t); |
94 | static void TF1Function(Int_t& unused1, Double_t* unused2, Double_t& chi2, Double_t *params, Int_t unused3); | |
e6339065 | 95 | |
96 | // static variable to be accessed by MINUIT | |
19442b86 | 97 | static TMatrixD* fgCorrelationMatrix; // contains fCurrentCorrelation in matrix form |
95e970ca | 98 | static TMatrixD* fgCorrelationMatrixSquared; // contains squared fCurrentCorrelation in matrix form |
19442b86 | 99 | static TMatrixD* fgCorrelationCovarianceMatrix; // contains the errors of fCurrentESD |
100 | static TVectorD* fgCurrentESDVector; // contains fCurrentESD | |
101 | static TVectorD* fgEntropyAPriori; // a-priori distribution for entropy regularization | |
95e970ca | 102 | static TVectorD* fgEfficiency; // efficiency |
9e065ad2 | 103 | /* |
95e970ca | 104 | static TVectorD* fgBinWidths; // bin widths to be taken into account in regularization |
9e065ad2 | 105 | static TVectorD* fgBinPos; // bin positions of unfolded |
106 | */ | |
107 | static TAxis *fgUnfoldedAxis; // bin widths and positions for unfolded | |
108 | static TAxis *fgMeasuredAxis; // bin widths and positions for measured | |
19442b86 | 109 | |
110 | static TF1* fgFitFunction; // fit function | |
e6339065 | 111 | |
19442b86 | 112 | // --- configuration params follow --- |
113 | static MethodType fgMethodType; // unfolding method to be used | |
114 | static Int_t fgMaxParams; // bins in unfolded histogram = number of fit params | |
115 | static Int_t fgMaxInput; // bins in measured histogram | |
116 | static Float_t fgOverflowBinLimit; // to fix fluctuations at high multiplicities, all entries above the limit are summarized in one bin | |
e6339065 | 117 | |
19442b86 | 118 | static RegularizationType fgRegularizationType; // regularization that is used during Chi2 method |
119 | static Float_t fgRegularizationWeight; // factor for regularization term | |
120 | static Int_t fgSkipBinsBegin; // (optional) skip the given number of bins in the regularization | |
121 | static Float_t fgMinuitStepSize; // (usually not needed to be changed) step size in minimization | |
9e065ad2 | 122 | static Float_t fgMinuitPrecision; // precision used by minuit. default = 1e-6 |
d22beb7f | 123 | static Int_t fgMinuitMaxIterations; // maximum number of iterations used by minuit. default = 5000 |
9f070ef5 | 124 | static Double_t fgMinuitStrategy; // minuit strategy: 0 (low), 1 (default), 2 (high) |
95e970ca | 125 | static Bool_t fgMinimumInitialValue; // set all initial values at least to the smallest value among the initial values |
126 | static Float_t fgMinimumInitialValueFix; // use this as the minimum initial value instead of determining it automatically | |
19442b86 | 127 | static Bool_t fgNormalizeInput; // normalize input spectrum |
95e970ca | 128 | static Float_t fgNotFoundEvents; // constraint on the total number of not found events sum(guess * (1/eff -1)) |
129 | static Bool_t fgSkipBin0InChi2; // skip bin 0 (= 0 measured) in chi2 function | |
e6339065 | 130 | |
19442b86 | 131 | static Float_t fgBayesianSmoothing; // smoothing parameter (0 = no smoothing) |
132 | static Int_t fgBayesianIterations; // number of iterations in Bayesian method | |
e6339065 | 133 | |
19442b86 | 134 | static Bool_t fgDebug; // debug flag |
135 | // --- end of configuration --- | |
95e970ca | 136 | |
137 | static Int_t fgCallCount; // call count to chi2 function | |
e6339065 | 138 | |
9e065ad2 | 139 | static Int_t fgPowern; // power of power law for regularization with power law |
140 | ||
141 | static Double_t fChi2FromFit; // Chi2 from fit at current iteration | |
142 | static Double_t fPenaltyVal; // Penalty value at current iteration (\beta * PU) | |
143 | static Double_t fAvgResidual; // Sum residuals / nbins | |
144 | ||
145 | static Int_t fgPrintChi2Details; // debug for chi2 calc | |
146 | ||
147 | // Pointers for interactive unfolder | |
148 | static TCanvas *fgCanvas; // Canvas for interactive unfolder | |
149 | static TH1 *fghUnfolded; // Unfolding result for interactive unfolder | |
150 | static TH2 *fghCorrelation; // Response matrix for interactive unfolder | |
151 | static TH1 *fghEfficiency; // Efficiency histo for interactive unfolder | |
152 | static TH1 *fghMeasured; // Measured distribution for interactive unfolder | |
153 | ||
19442b86 | 154 | private: |
e6339065 | 155 | AliUnfolding(const AliUnfolding&); |
156 | AliUnfolding& operator=(const AliUnfolding&); | |
157 | ||
158 | ClassDef(AliUnfolding, 0); | |
159 | }; | |
160 | ||
161 | #endif | |
162 |