Update to jet matching algorithm (from Redmer)
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskJetMatching.h
1 #ifndef AliAnalysisTaskJetMatching_H
2 #define AliAnalysisTaskJetMatching_H
3
4 #include <AliAnalysisTaskEmcalJetDev.h>
5 #include <AliEmcalJet.h>
6 #include <AliVTrack.h>
7 #include <TClonesArray.h>
8 #include <TMath.h>
9 #include <TRandom3.h>
10
11 class TF1;
12 class THF1;
13 class THF2;
14 class TProfile;
15 class AliRhoParameter;
16 class AliLocalRhoParameter;
17 class TClonesArray;
18
19 class AliAnalysisTaskJetMatching : public AliAnalysisTaskEmcalJetDev
20 {
21     public:
22         // enumerators
23         enum matchingSceme      {kGeoEtaPhi, kGeoR};
24         enum sourceBKG          {kNoSourceBKG, kSourceRho, kSourceLocalRho};
25         enum targetBKG          {kNoTargetBKG, kTargetRho, kTargetLocalRho};
26         // constructors, destructor
27                                 AliAnalysisTaskJetMatching();
28                                 AliAnalysisTaskJetMatching(const char *name);
29         virtual                 ~AliAnalysisTaskJetMatching();
30         // setting up the task and technical aspects
31         void                    ExecOnce();
32         virtual void            UserCreateOutputObjects();
33         TH1F*                   BookTH1F(const char* name, const char* x, Int_t bins, Double_t min, Double_t max, Bool_t append = kTRUE);
34         TH2F*                   BookTH2F(const char* name, const char* x, const char* y, Int_t binsx, Double_t minx, Double_t maxx, Int_t binsy, Double_t miny, Double_t maxy, Bool_t append = kTRUE);
35         virtual Bool_t          Run();
36         /* inline */    Double_t PhaseShift(Double_t x) const {  
37             while (x>=TMath::TwoPi())x-=TMath::TwoPi();
38             while (x<0.)x+=TMath::TwoPi();
39             return x; }
40         /* inline */    Double_t PhaseShift(Double_t x, Double_t n) const {
41             x = PhaseShift(x);
42             if(TMath::Nint(n)==2) while (x>TMath::Pi()) x-=TMath::Pi();
43             if(TMath::Nint(n)==3) {
44                 if(x>2.*TMath::TwoPi()/n) x = TMath::TwoPi() - x;
45                 if(x>TMath::TwoPi()/n) x = TMath::TwoPi()-(x+TMath::TwoPi()/n);
46             }
47             return x; }
48         /* inline */    Bool_t CompareTracks(AliVParticle* a, AliVParticle* b) const {
49             return (TMath::AreEqualAbs(a->Pt(), b->Pt(), 1e-2) && TMath::AreEqualAbs(a->Eta(), b->Eta(), 1e-2) && TMath::AreEqualAbs(a->Phi(), b->Phi(), 1e-2)) ? kTRUE : kFALSE; }
50         /* inline */    Double_t GetR(AliVParticle* a, AliVParticle* b) const {
51                if(!(a&&b)) return 999;
52                Double_t phiA(a->Phi()), phiB(b->Phi()), etaA(a->Eta()), etaB(b->Eta());
53                if(TMath::Abs(phiA-phiB) > TMath::Abs(phiA-phiB + TMath::TwoPi())) phiA+=TMath::TwoPi();
54                if(TMath::Abs(phiA-phiB) > TMath::Abs(phiA-phiB - TMath::TwoPi())) phiA-=TMath::TwoPi();
55                return TMath::Sqrt(TMath::Abs((etaA-etaB)*(etaA-etaB)+(phiA-phiB)*(phiA-phiB))); }
56
57         // setters - setup how to run
58         void                    SetDebugMode(Int_t d)                           {fDebug = d;}
59         void                    SetMatchingScheme(matchingSceme m)              {fMatchingScheme = m;}
60         void                    SetMatchConstituents(Bool_t m)                  {fMatchConstituents = m;}
61         void                    SetMinFracRecoveredConstituents(Float_t f)      {fMinFracRecoveredConstituents = f;}
62         void                    SetMinFracRecoveredConstituentPt(Float_t f)     {fMinFracRecoveredConstituentPt = f;}
63         void                    SetUseEmcalBaseJetCuts(Bool_t b)                {fUseEmcalBaseJetCuts = b;}
64         void                    SetGetBijection(Bool_t b)                       {fGetBijection = b;}
65         void                    SetSourceBKG(sourceBKG b)                       {fSourceBKG = b;}
66         void                    SetTargetBKG(targetBKG b)                       {fTargetBKG = b;}
67         void                    SetSourceJetsName(const char* n)                {fSourceJetsName = n;}
68         void                    SetTargetJetsName(const char* n)                {fTargetJetsName = n; }
69         void                    SetMatchedJetsName(const char* n)               {fMatchedJetsName = n;}
70         void                    SetSourceLocalRhoName(const char* n)            {fSourceRhoName = n;}
71         void                    SetTargetLocalRhoName(const char* n)            {fTargetRhoName = n;}
72         void                    SetSourceRadius(Float_t r)                      {fSourceRadius = r;}
73         void                    SetTargetRadius(Float_t r)                      {fTargetRadius = r;}
74         void                    SetMatchEta(Float_t f)                          {fMatchEta = f;}
75         void                    SetMatchPhi(Float_t f)                          {fMatchPhi = f;}
76         void                    SetMatchR(Float_t f)                            {fMatchR = f;}
77         // methods
78         void                    DoGeometricMatchingEtaPhi();
79         void                    DoGeometricMatchingR();
80         void                    DoConstituentMatching();
81         void                    GetBijection();
82         void                    PostMatchedJets();
83         // fill histogrmas
84         void                    FillAnalysisSummaryHistogram() const;
85         void                    FillMatchedJetHistograms() const;
86         // setters - analysis details
87         /* inline */    Bool_t PassesCuts(const AliVTrack* track) const {
88             return (!track || TMath::Abs(track->Eta()) > 0.9 || track->Phi() < 0 || track->Phi() > TMath::TwoPi()) ? kFALSE : kTRUE; }
89         /* inline */    Bool_t PassesCuts(AliEmcalJet* jet) const { return (jet) ? kTRUE : kFALSE; }
90         /* inline */    void ClearMatchedJetsCache() {
91             for (Int_t i(0); i < fNoMatchedJets; i++) {
92                 fMatchedJetContainer[i][0] = 0x0; fMatchedJetContainer[i][1] = 0x0; }
93             fNoMatchedJets = 0;
94         }
95         void                    PrintInfo() const;
96         virtual void            Terminate(Option_t* option);
97
98     private: 
99         Int_t                   fDebug;                 // debug level (0 none, 1 fcn calls, 2 verbose)
100         TClonesArray*           fSourceJets;            //! array with source jets
101         TString                 fSourceJetsName;        // name of array with source jets
102         TClonesArray*           fTargetJets;            //! array with target jets
103         TString                 fTargetJetsName;        // name of array with target jets
104         TClonesArray*           fMatchedJets;           //! final list of matched jets which is added to event
105         TString                 fMatchedJetsName;       // name of list of matched jets
106         AliLocalRhoParameter*   fSourceRho;             //! source rho
107         TString                 fSourceRhoName;         // source rho  name
108         AliLocalRhoParameter*   fTargetRho;             //! target rho
109         TString                 fTargetRhoName;         // target rho name
110         Bool_t                  fUseScaledRho;          // use scaled rho
111         Float_t                 fSourceRadius;          // source radius 
112         Float_t                 fTargetRadius;          // target radius
113         matchingSceme           fMatchingScheme;        // select your favorite matching algorithm
114         Bool_t                  fUseEmcalBaseJetCuts;   // use the emcal jet base class for jet cuts
115         sourceBKG               fSourceBKG;             // subtracted background of source jets
116         targetBKG               fTargetBKG;             // subtracted background of target jets
117         // transient object pointers
118         TList*                  fOutputList;            //! output list
119         TH1F*                   fHistUnsortedCorrelation;       //! dphi correlation of unsorted jets
120         TH1F*                   fHistMatchedCorrelation;        //! dphi correlation of matched jets
121         TH1F*                   fHistSourceJetPt;       //! pt of source jets
122         TH1F*                   fHistTargetJetPt;       //! pt of target jets
123         TH1F*                   fHistMatchedJetPt;      //! pt of matched jets
124         TH1F*                   fHistNoConstSourceJet;  //! number of constituents of source jet
125         TH1F*                   fHistNoConstTargetJet;  //! number of constituents of target jet
126         TH1F*                   fHistNoConstMatchJet;   //! number of constituents of matched jets
127         TProfile*               fProfFracPtMatched;     //! sum pt fraction for matched tracks / jet
128         TProfile*               fProfFracPtJets;        //! sum pt fraction for matched jets
129         TProfile*               fProfFracNoMatched;     //! no of constituents fraction found / jet
130         TProfile*               fProfFracNoJets;        //! no of consstituents fraction jet / jet
131         TH1F*                   fHistAnalysisSummary;   //! flags
132         TProfile*               fProfQAMatched;         //! QA spreads of matched jets
133         TProfile*               fProfQA;                //! QA spreads of source and target jets
134         Int_t                   fNoMatchedJets;         //! number of matched jets
135         AliEmcalJet*            fMatchedJetContainer[100][2];   //! pointers to matched jets
136         // geometric matching parameters
137         Float_t                 fMatchEta;              // max eta distance between centers of matched jets
138         Float_t                 fMatchPhi;              // max phi distance between centers of matched jets
139         Float_t                 fMatchR;                // max distance between matched jets
140         Bool_t                  fMatchConstituents;     // match constituents
141         Float_t                 fMinFracRecoveredConstituents;  // minimium fraction of constituents that needs to be found
142         Float_t                 fMinFracRecoveredConstituentPt; // minimium fraction of constituent pt that needs to be recovered
143         Bool_t                  fGetBijection;          // get bijection of source and matched jets
144
145         AliAnalysisTaskJetMatching(const AliAnalysisTaskJetMatching&);                  // not implemented
146         AliAnalysisTaskJetMatching& operator=(const AliAnalysisTaskJetMatching&);       // not implemented
147
148         ClassDef(AliAnalysisTaskJetMatching, 4);
149 };
150
151 #endif