26aa04d031557f3eef1555f10d4befa79fa7e252
[u/mrichter/AliRoot.git] / PWG / FLOW / Tasks / AliJetFlowTools.h
1 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * 
2  * See cxx source for full Copyright notice */ 
3  /* $Id$ */
4
5 /* Author: Redmer Alexander Bertens, rbertens@cern.ch, rbertens@nikhef.nl, r.a.bertens@uu.nl
6  * see implementation for additional information */
7
8 #ifndef AliJetFlowTools_H
9 #define AliJetFlowTools_H
10
11 // root forward declarations
12 class TF1;
13 class TF2;
14 class TH1D;
15 class TH2D;
16 class TCanvas;
17 class TString;
18 class TArrayD;
19 class TGraph;
20 class TGraphErrors;
21 class TObjArray;
22 // aliroot forward declarations
23 class AliAnaChargedJetResponseMaker;
24 class AliUnfolding;
25 // root includes
26 #include "TMatrixD.h"
27 #include "TList.h"
28 #include "TDirectoryFile.h"
29 #include "TFile.h"
30 #include "TProfile.h"
31 #include "TVirtualPad.h"
32 #include "TPaveText.h"
33 #include "TLegend.h"
34 #include "TLatex.h"
35 //_____________________________________________________________________________
36 class AliJetFlowTools {
37     public: 
38         AliJetFlowTools();
39     protected:
40         ~AliJetFlowTools();     // not implemented (deliberately). object ownership is a bit messy in this class
41                                 // since most (or all) of the objects are owned by the input and output files
42     public:
43         // enumerators
44         enum unfoldingAlgorithm {       // type of unfolding alrogithm
45             kChi2,                      // chi^2 unfolding, implemented in AliUnfolding
46             kBayesian,                  // Bayesian unfolding, implemented in RooUnfold
47             kBayesianAli,               // Bayesian unfolding, implemented in AliUnfolding
48             kSVD,                       // SVD unfolding, implemented in RooUnfold
49             kNone };                    // no unfolding
50         enum prior {                    // prior that is used for unfolding
51             kPriorChi2,                 // prior from chi^2 method
52             kPriorMeasured,             // use measured spectrum as prior
53             kPriorPythia };             // use pythia spectrum as prior
54         enum histoType {                // histogram identifier, only used internally
55             kInPlaneSpectrum,           // default style for spectrum
56             kOutPlaneSpectrum,
57             kUnfoldedSpectrum,
58             kFoldedSpectrum,
59             kMeasuredSpectrum,
60             kBar,                       // default style for bar histogram
61             kRatio,                     // default style for ratio
62             kV2,                        // default style for v2
63             kEmpty };                   // default style
64         // setters, interface to the class
65         void            SetSaveFull(Bool_t b)           {fSaveFull              = b;}
66         void            SetInputList(TList* list)       {
67             fInputList          = list;
68             fRefreshInput       = kTRUE;
69         }
70         void            SetOutputFileName(TString name) {fOutputFileName        = name;}
71         void            CreateOutputList(TString name) {
72             // create a new output list and add it to the full output
73             if(!fOutputFile) fOutputFile = new TFile(fOutputFileName.Data(), "RECREATE");
74             fOutputFile->cd();  // avoid nested dirs
75             if(name.EqualTo(fActiveString)) {
76                 printf(Form(" > Warning, duplicate output list, renaming %s to %s_2 ! < \n", name.Data(), name.Data()));
77                 name+="_2";
78             }
79             fActiveString = name;
80             fActiveDir = new TDirectoryFile(fActiveString.Data(), fActiveString.Data());
81             fActiveDir->cd();
82         }
83         void            SetCentralityBin(Int_t bin)             {
84             // in case of one centraltiy
85             fCentralityArray = new TArrayI(1);
86             fCentralityArray->AddAt(bin, 0);
87             // for one centrality there's no need for weights
88             fCentralityWeights = new TArrayD(1);
89             fCentralityWeights->AddAt(1., 0);
90         }
91         void            SetMergeSpectrumBins(TArrayI* a)        {fMergeBinsArray        =       a;}
92         void            SetCentralityBin(TArrayI* bins)         {
93             fCentralityArray = bins;
94         }
95         void            SetCentralityWeight(TArrayD* weights)   {
96             fCentralityWeights = weights;
97             if(!fCentralityArray) printf(" > Warning: centrality weights set, but bins are not defined! \n");
98         }
99         void            SetDetectorResponse(TH2D* dr)           {fDetectorResponse      = dr;}
100         void            SetJetFindingEfficiency(TH1D* e)        {fJetFindingEff         = e;}
101         void            SetBinsTrue(TArrayD* bins)              {fBinsTrue              = bins;}
102         void            SetBinsRec(TArrayD* bins)               {fBinsRec               = bins;}
103         void            SetBinsTruePrior(TArrayD* bins)         {fBinsTruePrior         = bins;}
104         void            SetBinsRecPrior(TArrayD* bins)          {fBinsRecPrior          = bins;}
105         void            SetSVDReg(Int_t r)                      {fSVDRegIn              = r; fSVDRegOut         = r;}
106         void            SetSVDReg(Int_t in, Int_t out)          {fSVDRegIn              = in; fSVDRegOut        = out;}
107         void            SetSVDToy(Bool_t b, Float_t r)          {fSVDToy                = b; fJetRadius         = r;}
108         void            SetBeta(Double_t b)                     {fBetaIn                = b; fBetaOut           = b;}
109         void            SetBeta(Double_t i, Double_t o)         {fBetaIn                = i; fBetaOut           = o;}
110         void            SetBayesianIter(Int_t i)                {fBayesianIterIn        = i; fBayesianIterOut    = i;}
111         void            SetBayesianIter(Int_t i, Int_t o)       {fBayesianIterIn        = i; fBayesianIterOut    = o;}
112         void            SetBayesianSmooth(Float_t s)            {fBayesianSmoothIn      = s; fBayesianSmoothOut = s;}
113         void            SetBayesianSmooth(Float_t i, Float_t o) {fBayesianSmoothIn      = i; fBayesianSmoothOut = o;}
114         void            SetAvoidRoundingError(Bool_t r)         {fAvoidRoundingError    = r;}
115         void            SetUnfoldingAlgorithm(unfoldingAlgorithm ua)    {fUnfoldingAlgorithm                    = ua;}
116         void            SetPrior(prior p)                       {fPrior                 = p;}
117         void            SetPrior(prior p, TH1D* spectrum)       {fPrior                 = p; fPriorUser         = spectrum;}
118         void            SetNormalizeSpectra(Bool_t b)           {fNormalizeSpectra      = b;}
119         void            SetNormalizeSpectra(Int_t e)            { // use to normalize to this no of events
120             fEventCount         = e;
121             fNormalizeSpectra   = kFALSE;
122         }
123         void            SetSmoothenPrior(Bool_t b, Float_t min = 50., Float_t max = 100., Float_t start= 75., Bool_t counts = kTRUE) {
124             fSmoothenPrior      = b;
125             fFitMin             = min;
126             fFitMax             = max;
127             fFitStart           = start;
128             fSmoothenCounts     = counts;
129         }
130         void            SetTestMode(Bool_t t)                   {fTestMode              = t;}
131         void            SetEventPlaneResolution(Double_t r)     {fEventPlaneRes         = r;}
132         void            SetUseDetectorResponse(Bool_t r)        {fUseDetectorResponse   = r;}
133         void            SetUseDptResponse(Bool_t r)             {fUseDptResponse        = r;}
134         void            SetTrainPowerFit(Bool_t t)              {fTrainPower            = t;}
135         void            SetDphiUnfolding(Bool_t i)              {fDphiUnfolding         = i;}
136         void            SetDphiDptUnfolding(Bool_t i)           {fDphiDptUnfolding      = i;}
137         void            SetExLJDpt(Bool_t i)                    {fExLJDpt               = i;}
138         void            SetWeightFunction(TF1* w)               {fResponseMaker->SetRMMergeWeightFunction(w);}
139         void            SetRMS(Bool_t r)                        {fRMS                   = r;}
140         void            SetSymmRMS(Bool_t r)                    {fSymmRMS               = r; fRMS               = r;}
141         void            Make();
142         void            MakeAU();       // test function, use with caution (09012014)
143         void            Finish() {
144             fOutputFile->cd();
145             if(fRMSSpectrumIn)  fRMSSpectrumIn->Write();
146             if(fRMSSpectrumOut) fRMSSpectrumOut->Write();
147             if(fRMSRatio)       fRMSRatio->Write();
148             fOutputFile->Close();}
149         void            PostProcess(
150                 TString def,
151                 Int_t columns = 4,
152                 Float_t rangeLow = 20,
153                 Float_t rangeUp = 80,
154                 TString in = "UnfoldedSpectra.root", 
155                 TString out = "ProcessedSpectra.root") const;
156         void            GetNominalValues(
157                 TH1D*& ratio,
158                 TGraphErrors*& v2,
159                 TArrayI* in,
160                 TArrayI* out,
161                 TString inFile = "UnfoldedSpectra.root",
162                 TString outFile = "Nominal.root") const;
163         void            GetCorrelatedUncertainty(
164                 TGraphAsymmErrors*& corrRatio,
165                 TGraphAsymmErrors*& corrV2,
166                 TArrayI* variationsIn,
167                 TArrayI* variationsOut,
168                 Bool_t sym,
169                 TArrayI* variantions2ndIn,
170                 TArrayI* variantions2ndOut,
171                 Bool_t sym2nd,
172                 TString type = "",
173                 TString type2 = "",
174                 Int_t columns = 4,
175                 Float_t rangeLow = 20,
176                 Float_t rangeUp = 80,
177                 Float_t corr = .5,
178                 TString in = "UnfoldedSpectra.root", 
179                 TString out = "CorrelatedUncertainty.root") const;
180         void            GetShapeUncertainty(
181                 TGraphAsymmErrors*& shapeRatio,
182                 TGraphAsymmErrors*& shapeV2,
183                 TArrayI* regularizationIn,
184                 TArrayI* regularizationOut,
185                 TArrayI* trueBinIn = 0x0,
186                 TArrayI* trueBinOut = 0x0,
187                 TArrayI* recBinIn = 0x0,
188                 TArrayI* recBinOut = 0x0,
189                 TArrayI* methodIn = 0x0,
190                 TArrayI* methodOut = 0x0,
191                 Int_t columns = 4,
192                 Float_t rangeLow = 20,
193                 Float_t rangeUp = 80,
194                 TString in = "UnfoldedSpectra.root", 
195                 TString out = "ShapeUncertainty.root") const;
196         Bool_t          SetRawInput (
197                 TH2D* detectorResponse, // detector response matrix
198                 TH1D* jetPtIn,          // in plane jet spectrum
199                 TH1D* jetPtOut,         // out of plane jet spectrum
200                 TH1D* dptIn,            // in plane delta pt distribution
201                 TH1D* dptOut,           // out of plane delta pt distribution
202                 Int_t eventCount = 0);  // event count (optional)
203         // static const helper functions, mainly histogram manipulation
204         static TH1D*    ResizeXaxisTH1D(TH1D* histo, Int_t low, Int_t up, TString suffix = "");
205         static TH2D*    ResizeYaxisTH2D(TH2D* histo, TArrayD* x, TArrayD* y, TString suffix = "");
206         static TH2D*    NormalizeTH2D(TH2D* histo, Bool_t noError = kTRUE);
207         static TH1D*    RebinTH1D(TH1D* histo, TArrayD* bins, TString suffix = "", Bool_t kill = kTRUE);
208         TH2D*           RebinTH2D(TH2D* histo, TArrayD* binsTrue, TArrayD* binsRec, TString suffix = "");
209         static TH2D*    MatrixMultiplication(TH2D* a, TH2D* b, TString name = "CombinedResponse");
210         static TH1D*    NormalizeTH1D(TH1D* histo, Double_t scale = 1.);
211         static TH1D*    MergeSpectrumBins(TArrayI* bins, TH1D* spectrum, TH2D* corr);
212         static TGraphErrors*    GetRatio(TH1 *h1 = 0x0, TH1* h2 = 0x0, TString name = "", Bool_t appendFit = kFALSE, Int_t xmax = -1);
213         static TGraphErrors*    GetV2(TH1* h1 = 0x0, TH1* h2 = 0x0, Double_t r = 0., TString name = "");
214         void     ReplaceBins(TArrayI* array, TGraphAsymmErrors* graph);
215         void     ReplaceBins(TArrayI* array, TGraphErrors* graph);
216         TGraphAsymmErrors*      GetV2WithSystematicErrors(
217                 TH1* h1, TH1* h2, Double_t r, TString name, 
218                 TH1* relativeErrorInUp,
219                 TH1* relativeErrorInLow,
220                 TH1* relativeErrorOutUp,
221                 TH1* relativeErrorOutLow,
222                 Float_t rho = 0.) const;
223         static void     GetSignificance(
224                 TGraphErrors* n,                // points with stat error
225                 TGraphAsymmErrors* shape,       // points with shape error
226                 TGraphAsymmErrors* corr,        // points with stat error
227                 Int_t low,                      // pt lower level
228                 Int_t up                        // pt upper level
229         );
230         static void     MinimizeChi22d();
231         static Double_t PhenixChi22d(const Double_t *xx );
232         static Double_t ConstructFunction2d(Double_t *x, Double_t *par);
233         static TF2*     ReturnFunction2d();
234         static void     MinimizeChi2nd();
235         static Double_t PhenixChi2nd(const Double_t *xx );
236         static Double_t ConstructFunctionnd(Double_t *x, Double_t *par);
237         static void     MinimizeChi2();
238         static TF2*     ReturnFunctionnd();
239         static Double_t PhenixChi2(const Double_t *xx );
240         static Double_t ConstructFunction(Double_t *x, Double_t *par);
241         static TF1*     ReturnFunction();
242         static void     WriteObject(TObject* object, TString suffix = "", Bool_t kill = kTRUE);
243         static TH2D*    ConstructDPtResponseFromTH1D(TH1D* dpt, Bool_t AvoidRoundingError);
244         static TH2D*    GetUnityResponse(TArrayD* binsTrue, TArrayD* binsRec, TString suffix = "");
245         void            SaveConfiguration(Bool_t convergedIn, Bool_t convergedOut) const;
246         static TMatrixD*        CalculatePearsonCoefficients(TMatrixD* covmat);
247         static TH1D*    SmoothenPrior(TH1D* spectrum, TF1* function, Double_t min, Double_t max, Double_t start, Bool_t kill = kTRUE, Bool_t counts = kTRUE);
248         // set style
249         void            SetTitleFontSize(Double_t s)    {fTitleFontSize = s;}
250         static void     Style(Bool_t legacy = kFALSE);
251         static void     Style(TCanvas* c, TString style = "PEARSON");
252         static void     Style(TVirtualPad* c, TString style = "SPECTRUM", Bool_t legacy = kFALSE);
253         static void     Style(TLegend* l);
254         static void     Style(TH1* h, EColor col = kBlue, histoType = kEmpty, Bool_t legacy = kFALSE);
255         static void     Style(TGraph* h, EColor col = kBlue, histoType = kEmpty, Bool_t legacy = kFALSE);
256         static TLegend* AddLegend(TVirtualPad* p, Bool_t style = kFALSE) {
257             if(!style) return p->BuildLegend(.565, .663, .882, .883);
258             else {
259                 TLegend* l = AddLegend(p, kFALSE);
260                 Style(l);
261                 return l;
262             }
263         }
264         static TPaveText*       AddTPaveText(
265                 // this text appears under the logo
266                 TString text, 
267                 Int_t r = 2,
268                 Double_t a = .587,
269                 Double_t b = .695,
270                 Double_t c = .872,
271                 Double_t d = .801) {
272             TPaveText* t(new TPaveText(a, b, c, d, "NDC"));
273             t->SetFillColor(0);            
274             t->SetBorderSize(0);
275             t->AddText(0.,0.,text.Data());
276             t->AddText(0., 0., Form("#it{R} = 0.%i anti-#it{k}_{T}, |#eta_{jet}|<%.1f", r, .9-r/10.));
277             t->SetTextColor(kBlack);
278             t->SetTextFont(42);
279             t->SetTextSize(gStyle->GetTextSize()*.8);
280             t->Draw("same");
281             return t;
282         } 
283         static TPaveText*       AddText(
284                 TString text, 
285                 EColor col,
286                 Double_t a = .2098,
287                 Double_t b = .5601,
288                 Double_t c = .613,
289                 Double_t d = .6211) {
290             TPaveText* t(new TPaveText(a, b, c, d, "NDC"));
291             t->SetFillColor(0);            
292             t->SetBorderSize(0);
293             t->AddText(0.,0.,text.Data());
294             t->SetTextColor(col);
295             t->SetTextFont(42);
296             t->SetTextSize(gStyle->GetTextSize()*.8);
297             t->Draw("same");
298             return t;
299         } 
300         static TLatex*          AddLogo(Int_t logo, Double_t xmin = .59, Double_t ymax = .81) {
301             if(logo == 0) return AddTLatex(xmin, ymax, "ALICE");
302             else if (logo == 1) return AddTLatex(xmin, ymax, "ALICE Preliminary");
303             else if (logo == 2) return AddTLatex(xmin, ymax, "ALICE Simulation");
304             return 0x0;
305         }
306         static TLatex*          AddSystem() {
307             return AddTLatex(0.55, 87, "Pb-Pb #sqrt{#it{s}}}_{NN} = 2.76 TeV");
308         }
309         static TLatex*          AddTLatex(Double_t xmin, Double_t ymax, TString string) {
310 TLatex* tex = new TLatex(xmin, ymax, string.Data());
311             tex->SetNDC();
312             tex->SetTextFont(42);
313             tex->Draw("same");
314             return tex;
315         }
316
317         static void     SavePadToPDF(TVirtualPad* pad)  {pad->SaveAs(Form("%s.pdf", pad->GetName()));}
318         // interface to AliUnfolding, not necessary but nice to have all parameters in one place
319         static void     SetMinuitStepSize(Float_t s)    {AliUnfolding::SetMinuitStepSize(s);}
320         static void     SetMinuitPrecision(Float_t s)   {AliUnfolding::SetMinuitPrecision(s);}
321         static void     SetMinuitPrecision(Int_t i)     {AliUnfolding::SetMinuitMaxIterations(i);}
322         static void     SetMinuitStrategy(Double_t s)   {AliUnfolding::SetMinuitStrategy(s);}
323         static void     SetDebug(Int_t d)               {AliUnfolding::SetDebug(d);}
324     private:
325         Bool_t          PrepareForUnfolding(); 
326         Bool_t          PrepareForUnfolding(Int_t low, Int_t up);
327         TH1D*           GetPrior(                       const TH1D* measuredJetSpectrum,
328                                                         const TH2D* resizedResponse,
329                                                         const TH1D* kinematicEfficiency,
330                                                         const TH1D* measuredJetSpectrumTrueBins,
331                                                         const TString suffix,
332                                                         const TH1D* jetFindingEfficiency);
333         TH1D*           UnfoldWrapper(                  const TH1D* measuredJetSpectrum, 
334                                                         const TH2D* resizedResponse,
335                                                         const TH1D* kinematicEfficiency,
336                                                         const TH1D* measuredJetSpectrumTrueBins,
337                                                         const TString suffix,
338                                                         const TH1D* jetFindingEfficiency = 0x0);
339         TH1D*           UnfoldSpectrumChi2(             const TH1D* measuredJetSpectrum, 
340                                                         const TH2D* resizedResponse,
341                                                         const TH1D* kinematicEfficiency,
342                                                         const TH1D* measuredJetSpectrumTrueBins,
343                                                         const TString suffix,
344                                                         const TH1D* jetFindingEfficiency = 0x0);
345         TH1D*           UnfoldSpectrumSVD(              const TH1D* measuredJetSpectrum, 
346                                                         const TH2D* resizedResponse,
347                                                         const TH1D* kinematicEfficiency,
348                                                         const TH1D* measuredJetSpectrumTrueBins,
349                                                         const TString suffix,
350                                                         const TH1D* jetFindingEfficiency = 0x0);
351         TH1D*           UnfoldSpectrumBayesianAli(      const TH1D* measuredJetSpectrum, 
352                                                         const TH2D* resizedResponse,
353                                                         const TH1D* kinematicEfficiency,
354                                                         const TH1D* measuredJetSpectrumTrueBins,
355                                                         const TString suffix,
356                                                         const TH1D* jetFindingEfficiency = 0x0);
357         TH1D*           UnfoldSpectrumBayesian(         const TH1D* measuredJetSpectrum, 
358                                                         const TH2D* resizedResponse,
359                                                         const TH1D* kinematicEfficiency,
360                                                         const TH1D* measuredJetSpectrumTrueBins,
361                                                         const TString suffix,
362                                                         const TH1D* jetFindingEfficiency = 0x0);
363         void            DoIntermediateSystematics(
364                 TArrayI* variationsIn,
365                 TArrayI* variationsOut,
366                 TH1D*& relativeErrorInUp,
367                 TH1D*& relativeErrorInLow,
368                 TH1D*& relativeErrorOutUp,
369                 TH1D*& relativeErrorOutLow,
370                 TH1D*& relativeSystematicIn,
371                 TH1D*& relativeSystematicOut,
372                 TH1D*& nominal,
373                 TH1D*& nominalIn,
374                 TH1D*& nominalOut,
375                 Int_t columns,
376                 Float_t rangeLow,
377                 Float_t rangeUp,
378                 TFile* readMe, 
379                 TString source = "",
380                 Bool_t RMS = kFALSE) const;
381         static void     ResetAliUnfolding();
382         // give object a unique name via the 'protect heap' functions. 
383         // may seem redundant, but some internal functions of root (e.g.
384         // ProjectionY()) check for existing objects by name and re-use them
385         TH1D*           ProtectHeap(TH1D* protect, Bool_t kill = kTRUE, TString suffix = "") const;
386         TH2D*           ProtectHeap(TH2D* protect, Bool_t kill = kTRUE, TString suffix = "") const;
387         TGraphErrors*   ProtectHeap(TGraphErrors* protect, Bool_t kill = kTRUE, TString suffix = "") const;
388         // members, accessible via setters
389         AliAnaChargedJetResponseMaker*  fResponseMaker; // utility object
390         Bool_t                  fRMS;                   // systematic method
391         Bool_t                  fSymmRMS;               // symmetric systematic method
392         TF1*                    fPower;                 // smoothening fit
393         Bool_t                  fSaveFull;              // save all generated histograms to file
394         TString                 fActiveString;          // identifier of active output
395         TDirectoryFile*         fActiveDir;             // active directory
396         TList*                  fInputList;             // input list
397         Bool_t                  fRefreshInput;          // re-read the input (called automatically if input list changes)
398         TString                 fOutputFileName;        // output file name
399         TFile*                  fOutputFile;            // output file
400         TArrayI*                fCentralityArray;       // array of centrality bins that are merged
401         TArrayI*                fMergeBinsArray;        // array of pt bins that are merged
402         TArrayD*                fCentralityWeights;     // array of centrality weights
403         TH2D*                   fDetectorResponse;      // detector response
404         TH1D*                   fJetFindingEff;         // jet finding efficiency
405         Double_t                fBetaIn;                // regularization strength, in plane unfolding
406         Double_t                fBetaOut;               // regularization strength, out of plane unfoldign
407         Int_t                   fBayesianIterIn;        // bayesian regularization parameter, in plane unfolding
408         Int_t                   fBayesianIterOut;       // bayesian regularization parameter, out plane unfolding
409         Float_t                 fBayesianSmoothIn;      // bayesian smoothening parameter (AliUnfolding)
410         Float_t                 fBayesianSmoothOut;     // bayesian smoothening parameter (AliUnfolding)
411         Bool_t                  fAvoidRoundingError;    // set dpt to zero for small values far from the diagonal
412         unfoldingAlgorithm      fUnfoldingAlgorithm;    // algorithm used for unfolding
413         prior                   fPrior;                 // prior for unfolding
414         TH1D*                   fPriorUser;             // user supplied prior (e.g. pythia spectrum)
415         TArrayD*                fBinsTrue;              // pt true bins
416         TArrayD*                fBinsRec;               // pt rec bins
417         TArrayD*                fBinsTruePrior;         // holds true bins for the chi2 prior for SVD. setting this is optional
418         TArrayD*                fBinsRecPrior;          // holds rec bins for the chi2 prior for SVD. setting this is optional
419         Int_t                   fSVDRegIn;              // svd regularization (a good starting point is half of the number of bins)
420         Int_t                   fSVDRegOut;             // svd regularization out of plane
421         Bool_t                  fSVDToy;                // use toy to estimate coveriance matrix for SVD method
422         Float_t                 fJetRadius;             // jet radius (for SVD toy)
423         Int_t                   fEventCount;            // number of events
424         Bool_t                  fNormalizeSpectra;      // normalize spectra to event count
425         Bool_t                  fSmoothenPrior;         // smoothen the tail of the measured spectrum using a powerlaw fit
426         Float_t                 fFitMin;                // lower bound of smoothening fit
427         Float_t                 fFitMax;                // upper bound of smoothening fit
428         Float_t                 fFitStart;              // from this value, use smoothening
429         Bool_t                  fSmoothenCounts;        // fill smoothened spectrum with counts
430         Bool_t                  fTestMode;              // unfold with unity response for testing
431         Bool_t                  fRawInputProvided;      // input histograms provided, not read from file
432         Double_t                fEventPlaneRes;         // event plane resolution for current centrality
433         Bool_t                  fUseDetectorResponse;   // add detector response to unfolding
434         Bool_t                  fUseDptResponse;        // add dpt response to unfolding
435         Bool_t                  fTrainPower;            // don't clear the params of fPower for call to Make
436                                                         // might give more stable results, but possibly introduces
437                                                         // a bias / dependency on previous iterations
438         Bool_t                  fDphiUnfolding;         // do the unfolding in in and out of plane orientation
439         Bool_t                  fDphiDptUnfolding;      // do the unfolding in dphi and dpt bins (to fit v2)
440         Bool_t                  fExLJDpt;               // exclude randon cones with leading jet
441         Double_t                fTitleFontSize;         // title font size
442         // members, set internally
443         TProfile*               fRMSSpectrumIn;         // rms of in plane spectra of converged unfoldings
444         TProfile*               fRMSSpectrumOut;        // rms of out of plane spectra of converged unfoldings
445         TProfile*               fRMSRatio;              // rms of ratio of converged unfolded spectra
446         TProfile*               fRMSV2;                 // rms of v2 of converged unfolded spectra
447         TH2D*                   fDeltaPtDeltaPhi;       // delta pt delta phi distribution
448         TH2D*                   fJetPtDeltaPhi;         // jet pt delta phi distribution
449         TH1D*                   fSpectrumIn;            // in plane jet pt spectrum
450         TH1D*                   fSpectrumOut;           // out of plane jet pt spectrum
451         TH1D*                   fDptInDist;             // in plane dpt distribution
452         TH1D*                   fDptOutDist;            // out of plane dpt distribution
453         TH2D*                   fDptIn;                 // in plane dpt matrix
454         TH2D*                   fDptOut;                // out plane dpt matrix
455         TH2D*                   fFullResponseIn;        // full response matrix, in plane
456         TH2D*                   fFullResponseOut;       // full response matrix, out of plane
457
458         // copy and assignment 
459         AliJetFlowTools(const AliJetFlowTools&);             // not implemented
460         AliJetFlowTools& operator=(const AliJetFlowTools&);  // not implemented
461 };
462 #endif
463 //_____________________________________________________________________________