PWGUD/base/AliMultiplicityCorrection.h

   1 /* $Id$ */
   2
   3 #ifndef ALIMULTIPLICITYCORRECTION_H
   4 #define ALIMULTIPLICITYCORRECTION_H
   5
   6 #include "TNamed.h"
   7
   8 //
   9 // class that contains the correction matrix and the functions for
  10 // correction the multiplicity spectrum
  11 // implements a several unfolding methods: e.g. chi2 minimization and bayesian unfolding
  12 //
  13
  14 class TH1;
  15 class TH2;
  16 class TH1F;
  17 class TH2F;
  18 class TH3F;
  19 class TF1;
  20 class TCollection;
  21
  22 // defined here, because it does not seem possible to predeclare these (or i do not know how)
  23 // -->
  24 // $ROOTSYS/include/TVectorDfwd.h:21: conflicting types for `typedef struct TVectorT<Double_t> TVectorD'
  25 // PWG0/dNdEta/AliMultiplicityCorrection.h:21: previous declaration as `struct TVectorD'
  26
  27 #include <TMatrixD.h>
  28 #include <TVectorD.h>
  29 #include <AliPWG0Helper.h>
  30 #include <AliUnfolding.h>
  31
  32 class AliMultiplicityCorrection : public TNamed {
  33   public:
  34     enum EventType { kTrVtx = 0, kMB, kINEL, kNSD };
  35     enum { kESDHists = 3, kMCHists = 4, kCorrHists = 6, kQualityRegions = 3 };
  36
  37     AliMultiplicityCorrection();
  38     AliMultiplicityCorrection(const Char_t* name, const Char_t* title);
  39     virtual ~AliMultiplicityCorrection();
  40
  41     static AliMultiplicityCorrection* Open(const char* fileName, const char* folderName = "Multiplicity");
  42
  43     virtual Long64_t Merge(const TCollection* list);
  44
  45     void FillMeasured(Float_t vtx, Int_t measured05, Int_t measured10, Int_t measured14);
  46     void FillTriggeredEvent(Int_t measured05, Int_t measured10, Int_t measured14);
  47     void FillGenerated(Float_t vtx, Bool_t triggered, Bool_t vertex, AliPWG0Helper::MCProcessType processType, Int_t generated05, Int_t generated10, Int_t generated14, Int_t generatedAll);
  48
  49     void FillCorrection(Float_t vtx, Int_t generated05, Int_t generated10, Int_t generated14, Int_t generatedAll, Int_t measured05, Int_t measured10, Int_t measured14);
  50     void FillNoVertexEvent(Float_t vtx, Bool_t vertexReconstructed, Int_t generated05, Int_t generated10, Int_t generated14, Int_t measured05, Int_t measured10, Int_t measured14);
  51
  52     Bool_t LoadHistograms(const Char_t* dir = 0);
  53     void SaveHistograms(const char* dir = 0);
  54     void DrawHistograms();
  55     void Rebin2DY(TH2F*& hist, Int_t nBins, Double_t* newBins) const;
  56     void Rebin3DY(TH3F*& hist, Int_t nBins, Double_t* newBins) const;
  57     void RebinGenerated(Int_t nBins05, Double_t* newBins05, Int_t nBins10, Double_t* newBins10, Int_t nBins13, Double_t* newBins13);
  58     void DrawComparison(const char* name, Int_t inputRange, Bool_t fullPhaseSpace, Bool_t normalizeESD, TH1* mcHist, Bool_t simple = kFALSE, EventType eventType = kTrVtx);
  59
  60     Int_t ApplyMinuitFit(Int_t inputRange, Bool_t fullPhaseSpace, EventType eventType, Int_t zeroBinEvents, Bool_t check = kFALSE, TH1* initialConditions = 0, Bool_t errorAsBias = kFALSE);
  61
  62     void ApplyBayesianMethod(Int_t inputRange, Bool_t fullPhaseSpace, EventType eventType, Float_t regPar = 1, Int_t nIterations = 100, TH1* initialConditions = 0, Int_t determineError = 1);
  63
  64     static TH1* CalculateStdDev(TH1** results, Int_t max);
  65     TH1* StatisticalUncertainty(AliUnfolding::MethodType methodType, Int_t inputRange, Bool_t fullPhaseSpace, EventType eventType, Int_t zeroBinEvents, Bool_t randomizeMeasured, Bool_t randomizeResponse, const TH1* compareTo = 0);
  66
  67     Int_t ApplyNBDFit(Int_t inputRange, Bool_t fullPhaseSpace, EventType eventType);
  68     void ApplyGaussianMethod(Int_t inputRange, Bool_t fullPhaseSpace);
  69
  70     void ApplyLaszloMethod(Int_t inputRange, Bool_t fullPhaseSpace, EventType eventType);
  71
  72     void Calculate0Bin(Int_t inputRange, EventType eventType, Int_t zeroBinEvents);
  73     Float_t GetFraction0Generated(Int_t inputRange);
  74
  75     TH2F* GetMultiplicityESD(Int_t i) const { return fMultiplicityESD[i]; }
  76     TH1F* GetTriggeredEvents(Int_t i) const { return fTriggeredEvents[i]; }
  77     TH2F* GetMultiplicityVtx(Int_t i) const { return fMultiplicityVtx[i]; }
  78     TH2F* GetMultiplicityMB(Int_t i) const { return fMultiplicityMB[i]; }
  79     TH2F* GetMultiplicityINEL(Int_t i) const { return fMultiplicityINEL[i]; }
  80     TH2F* GetMultiplicityNSD(Int_t i) const { return fMultiplicityNSD[i]; }
  81     TH2F* GetMultiplicityMC(Int_t i, EventType eventType) const;
  82     TH3F* GetCorrelation(Int_t i) const { return fCorrelation[i]; }
  83     TH1F* GetNoVertexEvents(Int_t i) const { return fNoVertexEvents[i]; }
  84     TH1F* GetMultiplicityESDCorrected(Int_t i) const { return fMultiplicityESDCorrected[i]; }
  85
  86     void SetMultiplicityESD(Int_t i, TH2F* const hist)  { fMultiplicityESD[i]  = hist; }
  87     void SetTriggeredEvents(Int_t i, TH1F* hist)  { fTriggeredEvents[i]  = hist; }
  88     void SetMultiplicityVtx(Int_t i, TH2F* const hist)  { fMultiplicityVtx[i]  = hist; }
  89     void SetMultiplicityMB(Int_t i, TH2F* const hist)   { fMultiplicityMB[i]   = hist; }
  90     void SetMultiplicityINEL(Int_t i, TH2F* const hist) { fMultiplicityINEL[i] = hist; }
  91     void SetMultiplicityNSD(Int_t i, TH2F* const hist) { fMultiplicityNSD[i] = hist; }
  92     void SetMultiplicityMC(Int_t i, EventType eventType, TH2F* const hist);
  93     void SetCorrelation(Int_t i, TH3F* const hist) { fCorrelation[i] = hist; }
  94     void SetNoVertexEvents(Int_t i, TH1F* hist) { fNoVertexEvents[i] = hist; }
  95     void SetMultiplicityESDCorrected(Int_t i, TH1F* const hist) { fMultiplicityESDCorrected[i] = hist; }
  96
  97     void SetGenMeasFromFunc(const TF1* inputMC, Int_t id);
  98     TH2F* CalculateMultiplicityESD(TH1* inputMC, Int_t correlationMap);
  99
 100     void FixTriggerEfficiencies(Int_t start);
 101
 102     TH1* GetConvoluted(Int_t i, EventType eventType);
 103     TH1* GetResiduals(Int_t i, EventType eventType, Float_t& residualSum);
 104
 105     void GetComparisonResults(Float_t* const mc = 0, Int_t* const mcLimit = 0, Float_t* const residuals = 0, Float_t* const ratioAverage = 0) const;
 106
 107     TH1* GetEfficiency(Int_t inputRange, EventType eventType);
 108     TH1* GetTriggerEfficiency(Int_t inputRange, EventType eventType);
 109
 110     static void SetQualityRegions(Bool_t SPDStudy);
 111     Float_t GetQuality(Int_t region) const { return fQuality[region]; }
 112
 113     void FFT(Int_t dir, Int_t m, Double_t *x, Double_t *y) const;
 114
 115     Float_t GetVertexBegin(Int_t i) { return fgVtxRangeBegin[i]; }
 116     Float_t GetVertexEnd(Int_t i) { return fgVtxRangeEnd[i]; }
 117
 118     void SetVertexRange(Int_t begin, Int_t end) { fVtxBegin = begin; fVtxEnd = end; }
 119
 120   protected:
 121     void SetupCurrentHists(Int_t inputRange, Bool_t fullPhaseSpace, EventType eventType);
 122
 123     Float_t BayesCovarianceDerivate(Float_t matrixM[251][251], const TH2* hResponse, Int_t k, Int_t i, Int_t r, Int_t u);
 124
 125     TH1* fCurrentESD;         //! current input esd
 126     TH2* fCurrentCorrelation; //! current correlation
 127     TH1* fCurrentEfficiency;  //! current efficiency
 128
 129     TH2F* fMultiplicityESD[kESDHists]; // multiplicity histogram: vtx vs multiplicity; array: |eta| < 0.5, 1.0, 1.4 (0..2)
 130     TH1F* fTriggeredEvents[kESDHists]; // (raw) multiplicity distribution of triggered events; array: |eta| < 0.5, 1.0, 1.4 (0..2)
 131     TH1F* fNoVertexEvents[kESDHists];  // distribution of true multiplicity just of triggered events without vertex or with 0 tracklets; array: |eta| < 0.5, 1.0, 1.4 (0..2)
 132
 133     TH2F* fMultiplicityVtx[kMCHists];  // multiplicity histogram of events that have a reconstructed vertex : vtx vs multiplicity; array: |eta| < 0.5, 1.0, 1.4, inf (0..3)
 134     TH2F* fMultiplicityMB[kMCHists];   // multiplicity histogram of triggered events                        : vtx vs multiplicity; array: |eta| < 0.5, 1.0, 1.4, inf (0..3)
 135     TH2F* fMultiplicityINEL[kMCHists]; // multiplicity histogram of all (inelastic) events                  : vtx vs multiplicity; array: |eta| < 0.5, 1.0, 1.4, inf (0..3)
 136     TH2F* fMultiplicityNSD[kMCHists]; // multiplicity histogram of NSD events                  : vtx vs multiplicity; array: |eta| < 0.5, 1.0, 1.4, inf (0..3)
 137
 138     TH3F* fCorrelation[kCorrHists];              // vtx vs. (gene multiplicity (trig+vtx)) vs. (meas multiplicity); array: |eta| < 0.5, 1, 1.4, (0..2 and 3..5), the first corrects to the eta range itself, the second to full phase space
 139
 140     TH1F* fMultiplicityESDCorrected[kCorrHists]; // corrected histograms
 141
 142     Int_t fLastBinLimit;        //! last bin limit, determined in SetupCurrentHists()
 143     Float_t fLastChi2MC;        //! last Chi2 between MC and unfolded ESD (calculated in DrawComparison)
 144     Int_t   fLastChi2MCLimit;   //! bin where the last chi2 breached a certain threshold, used to evaluate the multiplicity reach (calc. in DrawComparison)
 145     Float_t fLastChi2Residuals; //! last Chi2 of the ESD and the folded unfolded ESD (calculated in DrawComparison)
 146     Float_t fRatioAverage;      //! last average of |ratio-1| where ratio = unfolded / mc (bin 2..150)
 147
 148     Int_t fVtxBegin;            //! vertex range for analysis
 149     Int_t fVtxEnd;              //! vertex range for analysis
 150
 151     static Double_t fgVtxRangeBegin[kESDHists]; //! begin of allowed vertex range for this eta bin
 152     static Double_t fgVtxRangeEnd[kESDHists];   //! end of allowed vertex range for this eta bin
 153
 154     static Int_t   fgQualityRegionsB[kQualityRegions]; //! begin, given in multiplicity units
 155     static Int_t   fgQualityRegionsE[kQualityRegions]; //! end
 156     Float_t fQuality[kQualityRegions];                 //! stores the quality of the last comparison (calculated in DrawComparison). Contains 3 values that are averages of (MC - unfolded) / e(MC) in 3 regions, these are defined in fQualityRegionB,E
 157
 158  private:
 159     AliMultiplicityCorrection(const AliMultiplicityCorrection&);
 160     AliMultiplicityCorrection& operator=(const AliMultiplicityCorrection&);
 161
 162   ClassDef(AliMultiplicityCorrection, 7);
 163 };
 164
 165 #endif
 166