Line data Source code
1 : /* $Id$ */
2 :
3 : #ifndef ALIUNFOLDING_H
4 : #define ALIUNFOLDING_H
5 :
6 : //
7 : // class that implements several unfolding methods
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)
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;
23 : class TF1;
24 : class TCanvas;
25 : class TVirtualPad;
26 : class TAxis;
27 :
28 : class AliUnfolding : public TObject
29 : {
30 : public:
31 : enum RegularizationType { kNone = 0, kPol0, kPol1, kLog, kEntropy, kCurvature, kRatio, kPowerLaw, kLogLog };
32 : enum MethodType { kInvalid = -1, kChi2Minimization = 0, kBayesian = 1, kFunction = 2};
33 :
34 0 : virtual ~AliUnfolding() {};
35 :
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 0 : static void SetMinuitStepSize(Float_t stepSize) { fgMinuitStepSize = stepSize; }
42 0 : static void SetMinuitPrecision(Float_t pres) {fgMinuitPrecision = pres;}
43 0 : static void SetMinuitMaxIterations(Int_t iter) {fgMinuitMaxIterations = iter;}
44 0 : static void SetMinuitStrategy(Double_t strat) {fgMinuitStrategy = strat;}
45 0 : static void SetMinimumInitialValue(Bool_t flag, Float_t value = -1) { fgMinimumInitialValue = flag; fgMinimumInitialValueFix = value; }
46 0 : static void SetNormalizeInput(Bool_t flag) { fgNormalizeInput = flag; }
47 0 : static void SetNotFoundEvents(Float_t notFoundEvents) { fgNotFoundEvents = notFoundEvents; }
48 0 : static void SetSkip0BinInChi2(Bool_t flag) { fgSkipBin0InChi2 = flag; }
49 : static void SetBayesianParameters(Float_t smoothing, Int_t nIterations);
50 : static void SetFunction(TF1* function);
51 0 : static void SetDebug(Bool_t flag) { fgDebug = flag; }
52 0 : static void SetPowern(Int_t n) {fgPowern = -1*n;}
53 :
54 : static Int_t Unfold(TH2* correlation, TH1* efficiency, TH1* measured, TH1* initialConditions, TH1* result, Bool_t check = kFALSE);
55 : static Int_t UnfoldGetBias(TH2* correlation, TH1* efficiency, TH1* measured, TH1* initialConditions, TH1* result);
56 :
57 : static TH1* GetPenaltyPlot(Double_t* params);
58 : static TH1* GetPenaltyPlot(TH1* corrected);
59 :
60 : static TH1* GetResidualsPlot(Double_t* params);
61 : static TH1* GetResidualsPlot(TH1* corrected);
62 :
63 0 : static Double_t GetChi2FromFit() {return fChi2FromFit;}
64 0 : static Double_t GetPenaltyVal() {return fPenaltyVal;}
65 0 : 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 :
71 : protected:
72 : AliUnfolding() {};
73 :
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 :
78 : static void CreateOverflowBin(TH2* correlation, TH1* measured);
79 : static void SetStaticVariables(TH2* correlation, TH1* measured, TH1* efficiency);
80 :
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 :
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);
89 : static Double_t RegularizationRatio(TVectorD& params);
90 : static Double_t RegularizationPowerLaw(TVectorD& params);
91 : static Double_t RegularizationLogLog(TVectorD& params);
92 :
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);
95 :
96 : // static variable to be accessed by MINUIT
97 : static TMatrixD* fgCorrelationMatrix; // contains fCurrentCorrelation in matrix form
98 : static TMatrixD* fgCorrelationMatrixSquared; // contains squared fCurrentCorrelation in matrix form
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
102 : static TVectorD* fgEfficiency; // efficiency
103 : /*
104 : static TVectorD* fgBinWidths; // bin widths to be taken into account in regularization
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
109 :
110 : static TF1* fgFitFunction; // fit function
111 :
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
117 :
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
122 : static Float_t fgMinuitPrecision; // precision used by minuit. default = 1e-6
123 : static Int_t fgMinuitMaxIterations; // maximum number of iterations used by minuit. default = 5000
124 : static Double_t fgMinuitStrategy; // minuit strategy: 0 (low), 1 (default), 2 (high)
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
127 : static Bool_t fgNormalizeInput; // normalize input spectrum
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
130 :
131 : static Float_t fgBayesianSmoothing; // smoothing parameter (0 = no smoothing)
132 : static Int_t fgBayesianIterations; // number of iterations in Bayesian method
133 :
134 : static Bool_t fgDebug; // debug flag
135 : // --- end of configuration ---
136 :
137 : static Int_t fgCallCount; // call count to chi2 function
138 :
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 :
154 : private:
155 : AliUnfolding(const AliUnfolding&);
156 : AliUnfolding& operator=(const AliUnfolding&);
157 :
158 170 : ClassDef(AliUnfolding, 0);
159 : };
160 :
161 : #endif
162 :
|
