]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskJetMatching.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskJetMatching.cxx
CommitLineData
7f7120c5 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
f1354962 16//
17// General task to match two sets of jets
18//
19// This task takes two TClonesArray's as input:
20// [1] fSourceJets - e.g. pythia jets
21// [2] fTargetJets - e.g. a samle containing pythia jets embedded in a pbpb event
22// the task will try to match jets from the source array to the target array, and
23// save the found TARGET jets in a new array 'fMatchedJets'
c2460f8a 24// the purpose of this task is twofold
25// [1] matching for embedding studies
26// [2] matching to create a detector response function
f1354962 27//
285db795 28// matching is done in three steps
29// [1] geometric matching, where jets are matched by R = sqrt(dphi*dphi-deta*deta) or directly via dphi and deta
30// [2] optional injection / bijection check
31// in this check, it is made sure that fSourceJets -> fMatchedJets is either an injective non-surjective
32// or bijective map, depending on the matching resolution which is chosen for step [1]
33// so that each source jet has a unique target and vice-versa.
34// if R (or dphi, deta) are proportional to, or larger than, the jet radius, matching tends to be
35// bijective (each source has a target), if R is chosen to be very small, source jets might be lost
36// as no target can be found.
37// the mapping is done in such a way that each target is matched to its closest source and each source
38// is mapped to its closest target
39// [3] constituent matching
40// - how many constituents of the source jet are present in the matched jet?
41// a cut on the constituent fraction can be performed (not recommended)
42// - how much of the original pt is recovered in the matched jet?
43// a cut on the fraction of recovered / original pt can be performed (recommended)
44//
c2460f8a 45// detector response
46// to obtain a detector respose function, supply
47// [1] fSourceJets - particle level jets
48// [2] fTargetJets - detector level jets
49//
f1354962 50// Author: Redmer Alexander Bertens, Utrecht Univeristy, Utrecht, Netherlands
51// rbertens@cern.ch, rbertens@nikhef.nl, r.a.bertens@uu.nl
52
53// root includes
54#include <TClonesArray.h>
55#include <TChain.h>
56#include <TMath.h>
57#include <TF1.h>
58#include <TF2.h>
10d7808d 59#include <TH3.h>
f1354962 60#include <TH1F.h>
61#include <TH2F.h>
10d7808d 62#include <TH3F.h>
f1354962 63#include <TProfile.h>
7ec64e26 64#include <TFile.h>
65#include <TTree.h>
66#include <TKey.h>
67#include <TSystem.h>
f1354962 68// aliroot includes
69#include <AliAnalysisTask.h>
70#include <AliAnalysisManager.h>
71#include <AliLog.h>
72#include <AliVEvent.h>
73#include <AliVParticle.h>
74// emcal jet framework includes
75#include <AliEmcalJet.h>
76#include <AliAnalysisTaskJetMatching.h>
f6d1b1a7 77#include <AliLocalRhoParameter.h>
f1354962 78
79class AliAnalysisTaskJetMatching;
80using namespace std;
81
82ClassImp(AliAnalysisTaskJetMatching)
83
9239b066 84AliAnalysisTaskJetMatching::AliAnalysisTaskJetMatching() : AliAnalysisTaskEmcalJet("AliAnalysisTaskJetMatching", kTRUE),
90c6bce5 85 fSourceJets(0), fSourceJetsName(0), fTargetJets(0), fTargetJetsName(0), fMatchedJets(0), fMatchedJetsName(GetName()), fSourceRho(0), fSourceRhoName(0), fTargetRho(0), fTargetRhoName(0), fUseScaledRho(0), fSourceRadius(0.3), fTargetRadius(0.3), fMatchingScheme(kGeoEtaPhi), fUseEmcalBaseJetCuts(kFALSE), fSourceBKG(kNoSourceBKG), fTargetBKG(kNoTargetBKG), fOutputList(0), fHistUnsortedCorrelation(0), fHistMatchedCorrelation(0), fHistSourceJetPt(0), fHistMatchedSourceJetPt(0), fHistTargetJetPt(0), fHistMatchedJetPt(0), fHistSourceMatchedJetPt(0), fHistDetectorResponseProb(0), fHistNoConstSourceJet(0), fHistNoConstTargetJet(0), fHistNoConstMatchJet(0), fProfFracPtMatched(0), fProfFracPtJets(0), fProfFracNoMatched(0), fProfFracNoJets(0), fHistDiJet(0), fHistDiJetLeadingJet(0), fHistDiJetDPhi(0), fHistDiJetDPt(0), fHistAnalysisSummary(0), fProfQAMatched(0), fProfQA(0), fNoMatchedJets(200), fMatchEta(.3), fMatchPhi(.3), fMatchR(.08), fDoDetectorResponse(kFALSE), fMatchConstituents(kTRUE), fMinFracRecoveredConstituents(.5), fMinFracRecoveredConstituentPt(0.5), fGetBijection(kTRUE), fh1Trials(0x0), fh1AvgTrials(0x0), fh1Xsec(0x0), fAvgTrials(0) {
f1354962 86 // default constructor
87 ClearMatchedJetsCache();
88}
89//_____________________________________________________________________________
9239b066 90AliAnalysisTaskJetMatching::AliAnalysisTaskJetMatching(const char* name) : AliAnalysisTaskEmcalJet(name, kTRUE),
90c6bce5 91 fSourceJets(0), fSourceJetsName(0), fTargetJets(0), fTargetJetsName(0), fMatchedJets(0), fMatchedJetsName(GetName()), fSourceRho(0), fSourceRhoName(0), fTargetRho(0), fTargetRhoName(0), fUseScaledRho(0), fSourceRadius(0.3), fTargetRadius(0.3), fMatchingScheme(kGeoEtaPhi), fUseEmcalBaseJetCuts(kFALSE), fSourceBKG(kNoSourceBKG), fTargetBKG(kNoTargetBKG), fOutputList(0), fHistUnsortedCorrelation(0), fHistMatchedCorrelation(0), fHistSourceJetPt(0), fHistMatchedSourceJetPt(0), fHistTargetJetPt(0), fHistMatchedJetPt(0), fHistSourceMatchedJetPt(0), fHistDetectorResponseProb(0), fHistNoConstSourceJet(0), fHistNoConstTargetJet(0), fHistNoConstMatchJet(0), fProfFracPtMatched(0), fProfFracPtJets(0), fProfFracNoMatched(0), fProfFracNoJets(0), fHistDiJet(0), fHistDiJetLeadingJet(0), fHistDiJetDPhi(0), fHistDiJetDPt(0), fHistAnalysisSummary(0), fProfQAMatched(0), fProfQA(0), fNoMatchedJets(200), fMatchEta(.3), fMatchPhi(.3), fMatchR(.08), fDoDetectorResponse(kFALSE), fMatchConstituents(kTRUE), fMinFracRecoveredConstituents(0.5), fMinFracRecoveredConstituentPt(0.5), fGetBijection(kTRUE), fh1Trials(0x0), fh1AvgTrials(0x0), fh1Xsec(0x0), fAvgTrials(0) {
f1354962 92 // constructor
93 ClearMatchedJetsCache();
94 DefineInput(0, TChain::Class());
95 DefineOutput(1, TList::Class());
96}
97//_____________________________________________________________________________
98AliAnalysisTaskJetMatching::~AliAnalysisTaskJetMatching()
99{
100 // destructor
101 if(fOutputList) delete fOutputList;
102}
103//_____________________________________________________________________________
104void AliAnalysisTaskJetMatching::ExecOnce()
105{
106 // initialize the anaysis
90c6bce5 107 #ifdef DEBUGTASK
108 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
109 #endif
f1354962 110 // get the stand alone jets from the input event
111 fSourceJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fSourceJetsName.Data()));
112 if(!fSourceJets) AliFatal(Form("%s: Container with name %s not found. Aborting", GetName(), fSourceJetsName.Data()));
113 fTargetJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTargetJetsName.Data()));
114 if(!fTargetJets) AliFatal(Form("%s: Container with name %s not found. Aborting", GetName(), fSourceJetsName.Data()));
115 // append the list of matched jets to the event
116 fMatchedJets->Delete();
117 if (!(InputEvent()->FindListObject(fMatchedJetsName))) InputEvent()->AddObject(fMatchedJets);
118 else AliFatal(Form("%s: Object with name %s already in event! Aborting", GetName(), fMatchedJetsName.Data()));
119 FillAnalysisSummaryHistogram();
f6d1b1a7 120 switch (fSourceBKG) {
121 case kSourceLocalRho : {
122 fSourceRho = dynamic_cast<AliLocalRhoParameter*>(InputEvent()->FindListObject(fSourceRhoName));
123 if(!fSourceRho) AliFatal(Form("%s: Object with name %s requested but not found! Aborting", GetName(), fSourceRhoName.Data()));
124 } break;
125 default : break;
126 }
127 switch (fTargetBKG) {
128 case kTargetLocalRho : {
129 fTargetRho = dynamic_cast<AliLocalRhoParameter*>(InputEvent()->FindListObject(fTargetRhoName));
130 if(!fTargetRho) AliFatal(Form("%s: Object with name %s requested but not found! Aborting", GetName(), fTargetRhoName.Data()));
131 } break;
132 default : break;
133 }
9239b066 134 AliAnalysisTaskEmcalJet::ExecOnce(); // init base class
c2460f8a 135 if(fDoDetectorResponse) fMatchConstituents = kFALSE;
f1354962 136}
137//_____________________________________________________________________________
138void AliAnalysisTaskJetMatching::UserCreateOutputObjects()
139{
140 // create output objects
90c6bce5 141 #ifdef DEBUGTASK
142 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
143 #endif
144 // add the matched jets array to the event
f1354962 145 fMatchedJets = new TClonesArray("AliEmcalJet");
146 fMatchedJets->SetName(fMatchedJetsName);
147 fOutputList = new TList();
148 fOutputList->SetOwner(kTRUE);
149 // add analysis histograms
150 fHistUnsortedCorrelation = BookTH1F("fHistUnsortedCorrelation", "#Delta #varphi unsorted", 50, 0, TMath::Pi());
151 fHistMatchedCorrelation = BookTH1F("fHistMatchedCorrelation", "#Delta #varphi matched", 50, 0, TMath::Pi());
2e3ee924 152 fHistSourceJetPt = (fDoDetectorResponse) ? BookTH1F("fHistParticleLevelJetPt", "p_{t}^{gen} [GeV/c]", 150, -50, 250) : BookTH1F("fHistSourceJetPt", "p_{t} [GeV/c]", 150, -50, 250);
153 fHistMatchedSourceJetPt = (fDoDetectorResponse) ? BookTH1F("fHistMatchedParticleLevelJetPt", "p_{t}^{gen} [GeV/c]", 150, -50, 250) : BookTH1F("fHistMatchedSourceJetPt", "p_{t} [GeV/c]", 150, -50, 250);
c2460f8a 154 fHistTargetJetPt = BookTH1F("fHistTargetJetPt", "p_{t} [GeV/c]", 150, -50, 250);
2e3ee924 155 fHistMatchedJetPt = (fDoDetectorResponse) ? BookTH1F("fHistDetectorLevelJet", "p_{t}^{rec}", 150, -50, 250) : BookTH1F("fHistMatchedJetPt", "p_{t} [GeV/c]", 150, -50, 250);
7f7120c5 156 fHistSourceMatchedJetPt = (fDoDetectorResponse) ? BookTH2F("fHistDetectorResponse", "particle level jet p_{t}^{gen} [GeV/c]", "detector level jet p_{t}^{rec} [GeV/c]", 300, -50, 250, 300, -50, 250) : BookTH2F("fHistSourceMatchedJetPt", "source jet p_{t} [GeV/c]", "matched jet p_{t} [GeV/c]", 300, -50, 250, 300, -50, 250);
cd5876d1 157 if(fDoDetectorResponse) {
7ec64e26 158 fHistDetectorResponseProb = BookTH2F("fHistDetectorResponseProb", "(p_{t}^{det} - p_{t}^{part})/p_{t}^{part}", "p_{t}^{part}", 100, -1.5, 1., 20, 0, 200);
159 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
160 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
161 fOutputList->Add(fh1Xsec);
162 fh1Trials = new TH1F("fh1Trials","trials root file",1,0,1);
163 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
164 fOutputList->Add(fh1Trials);
165 fh1AvgTrials = new TH1F("fh1AvgTrials","avg trials root file",1,0,1);
166 fh1AvgTrials->GetXaxis()->SetBinLabel(1,"#sum{avg ntrials}");
167 fOutputList->Add(fh1AvgTrials);
cd5876d1 168 }
7f7120c5 169 fHistNoConstSourceJet = BookTH1F("fHistNoConstSourceJet", "number of constituents source jets", 50, 0, 50);
170 fHistNoConstTargetJet = BookTH1F("fHistNoConstTargetJet", "number of constituents target jets", 50, 0, 50);
171 fHistNoConstMatchJet = BookTH1F("fHistNoConstMatchJet", "number of constituents matched jets", 50, 0, 50);
2e3ee924 172 if(!fDoDetectorResponse) { // these observables cannot be measured in current detector response mode
7f7120c5 173 fProfFracPtMatched = new TProfile("fProfFracPtMatched", "recovered target p_{T} / source p_{T}", 15, -50, 250);
2e3ee924 174 fOutputList->Add(fProfFracPtMatched);
7f7120c5 175 fProfFracNoMatched = new TProfile("fProfFracNoMatched", "recovered target constituents / source constituents", 15, -50, 250);
2e3ee924 176 fOutputList->Add(fProfFracNoMatched);
177 }
f1354962 178 // the analysis summary histo which stores all the analysis flags is always written to file
179 fHistAnalysisSummary = BookTH1F("fHistAnalysisSummary", "flag", 51, -0.5, 15.5);
180 fProfQAMatched = new TProfile("fProfQAMatched", "fProfQAMatched", 3, 0, 3);
181 fProfQAMatched->GetXaxis()->SetBinLabel(1, "<#delta p_{t} >");
182 fProfQAMatched->GetXaxis()->SetBinLabel(2, "<#delta #eta>");
183 fProfQAMatched->GetXaxis()->SetBinLabel(3, "<#delta #varphi>");
184 fOutputList->Add(fProfQAMatched);
185 fProfQA = new TProfile("fProfQA", "fProfQA", 3, 0, 3);
186 fProfQA->GetXaxis()->SetBinLabel(1, "<#delta p_{t} >");
187 fProfQA->GetXaxis()->SetBinLabel(2, "<#delta #eta>");
188 fProfQA->GetXaxis()->SetBinLabel(3, "<#delta #varphi>");
189 fOutputList->Add(fProfQA);
7f7120c5 190 fProfFracPtJets = new TProfile("fProfFracPtJets", "source p_{T} / target p_{T}", 15, -50, 250);
0aaf31ab 191 fOutputList->Add(fProfFracPtJets);
7f7120c5 192 fProfFracNoJets = new TProfile("fProfFracNoJets", "source constituents / target constituents", 15, -50, 250);
0aaf31ab 193 fOutputList->Add(fProfFracNoJets);
08fae484 194 switch (fMatchingScheme) {
195 case kDiJet : {
10d7808d 196 fHistDiJet = BookTH3F("fHistDiJet", "matched di-jet #varphi", "matched di-jet #eta", "leading jet p_{T} (GeV/c)", 100, 0., TMath::TwoPi(), 100, -5., 5., 100, 0, 200);
197 fHistDiJetLeadingJet = BookTH3F("fHistDiJetLeadingJet", "leading jet #varphi", "leadingd jet #eta", "leading jet p_{T} (GeV/c)", 100, 0., TMath::TwoPi(), 100, -5., 5., 100, 0, 200);
198 fHistDiJetDPhi = BookTH2F("fHistDiJetDPhi", "leading jet #varphi - (matched jet #varphi - #pi)", "leading jet p_{T} (GeV/c)", 100, -1.*TMath::Pi(), TMath::Pi(), 100, 0, 200);
199 fHistDiJetDPt = BookTH2F("fHistDiJetDPt", "leading jet p_{T} - sub leading jet p_{T} (GeV/c)", "leading jet p_{T} (GeV/c)", 100, -25, 25, 100, 0, 200);
08fae484 200 } break;
201 default : break;
202 }
f1354962 203 fOutputList->Sort();
204 PostData(1, fOutputList);
205}
206//_____________________________________________________________________________
207TH1F* AliAnalysisTaskJetMatching::BookTH1F(const char* name, const char* x, Int_t bins, Double_t min, Double_t max, Bool_t append)
208{
209 // book a TH1F and connect it to the output container
90c6bce5 210 #ifdef DEBUGTASK
211 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
212 #endif
10d7808d 213 if(append && !fOutputList) return 0x0;
f1354962 214 TString title(name);
215 title += Form(";%s;[counts]", x);
216 TH1F* histogram = new TH1F(name, title.Data(), bins, min, max);
217 histogram->Sumw2();
218 if(append) fOutputList->Add(histogram);
219 return histogram;
220}
221//_____________________________________________________________________________
222TH2F* AliAnalysisTaskJetMatching::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)
223{
224 // book a TH2F and connect it to the output container
90c6bce5 225 #ifdef DEBUGTASK
226 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
227 #endif
10d7808d 228 if(append && !fOutputList) return 0x0;
f1354962 229 TString title(name);
230 title += Form(";%s;%s", x, y);
231 TH2F* histogram = new TH2F(name, title.Data(), binsx, minx, maxx, binsy, miny, maxy);
232 histogram->Sumw2();
233 if(append) fOutputList->Add(histogram);
234 return histogram;
235}
236//_____________________________________________________________________________
10d7808d 237TH3F* AliAnalysisTaskJetMatching::BookTH3F(const char* name, const char* x, const char* y, const char* z, Int_t binsx, Double_t minx, Double_t maxx, Int_t binsy, Double_t miny, Double_t maxy, Int_t binsz, Double_t minz, Double_t maxz, Bool_t append)
238{
239 // book a TH2F and connect it to the output container
90c6bce5 240 #ifdef DEBUGTASK
241 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
242 #endif
10d7808d 243 if(append && !fOutputList) return 0x0;
244 TString title(name);
245 title += Form(";%s;%s;%s", x, y, z);
246 TH3F* histogram = new TH3F(name, title.Data(), binsx, minx, maxx, binsy, miny, maxy, binsz, minz, maxz);
247 histogram->Sumw2();
248 if(append) fOutputList->Add(histogram);
249 return histogram;
250}
251//_____________________________________________________________________________
7ec64e26 252Bool_t AliAnalysisTaskJetMatching::Notify()
253{
254 // for each file get the number of trials and pythia cross section
255 // see $ALICE_ROOT/PWGJE/AliAnalysisTaskJetProperties.cxx
256 // $ALICE_ROOT/PWG/Tools/AliAnalysisHelperJetTasks.cxx
257 // this function is implenented here temporarily to avoid introducing a dependency
258 // later on this could just be a call to a static helper function
90c6bce5 259 #ifdef DEBUGTASK
260 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
261 #endif
7ec64e26 262 Float_t xsection(0), ftrials(1);
263 fAvgTrials = ftrials;
264 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
265 if(tree) {
266 TFile *curfile = tree->GetCurrentFile();
267 if (!curfile) return kFALSE;
268 TString file(curfile->GetName());
269 if(file.Contains("root_archive.zip#")) file.Replace(file.Index("#",1,file.Index("root_archive",12,0,TString::kExact),TString::kExact)+1,file.Index(".root",5,TString::kExact)-file.Index("root_archive",12,0,TString::kExact),"");
270 else file.ReplaceAll(gSystem->BaseName(file.Data()),"");
271 TFile *fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec.root"));
272 if(!fxsec) {
273 fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
274 if(fxsec) {
275 TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
276 if(key) {
277 TList *list = dynamic_cast<TList*>(key->ReadObj());
278 if(list) {
279 xsection = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
280 ftrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
281 }
282 }
283 fxsec->Close();
284 }
285 } else {
286 TTree *xtree = (TTree*)fxsec->Get("Xsection");
287 if(xtree) {
288 UInt_t _ftrials = 0;
289 Double_t _xsection = 0;
290 xtree->SetBranchAddress("xsection",&_xsection);
291 xtree->SetBranchAddress("ntrials",&_ftrials);
292 xtree->GetEntry(0);
293 ftrials = _ftrials;
294 xsection = _xsection;
295 }
296 fxsec->Close();
297 }
298 fh1Xsec->Fill("<#sigma>",xsection);
299 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
300 if(ftrials >= nEntries && nEntries > 0.) fAvgTrials = ftrials/nEntries;
301 fh1Trials->Fill("#sum{ntrials}",ftrials);
302 }
303 return kTRUE;
304}
305//_____________________________________________________________________________
f1354962 306Bool_t AliAnalysisTaskJetMatching::Run()
307{
308 // execute once for each event
90c6bce5 309 #ifdef DEBUGTASK
310 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
311 #endif
459e1638 312 if(!(InputEvent() && fSourceJets && fTargetJets && IsEventSelected())) return kFALSE;
7ec64e26 313 if(fh1AvgTrials) fh1AvgTrials->Fill("#sum{avg ntrials}", fAvgTrials);
285db795 314 // step one: do geometric matching
f1354962 315 switch (fMatchingScheme) {
08fae484 316 // FIXME having separate dedicated functions is not necessary and historical
317 // cluttered code - should be cleaned up and merged into one
f1354962 318 case kGeoEtaPhi : {
319 DoGeometricMatchingEtaPhi();
c52fd114 320 // break if no jet was matched
321 if(!fMatchedJetContainer[1][0]) return kTRUE;
f1354962 322 } break;
323 case kGeoR : {
324 DoGeometricMatchingR();
c52fd114 325 // break if no jet was matched
326 if(!fMatchedJetContainer[1][0]) return kTRUE;
f1354962 327 } break;
08fae484 328 case kDiJet : {
329 DoDiJetMatching();
330 } break;
285db795 331 default : break;
f1354962 332 }
285db795 333 // optional step two: get a bijection (avoid duplicate matches)
334 if(fGetBijection) GetBijection();
335 // optional step three: match constituents within matched jets
336 if(fMatchConstituents) DoConstituentMatching();
f1354962 337 // stream data to output
338 PostMatchedJets();
339 FillMatchedJetHistograms();
90c6bce5 340 #ifdef DEBUGTASK
341 PrintInfo();
342 #endif
f1354962 343 PostData(1, fOutputList);
344 return kTRUE;
345}
346//_____________________________________________________________________________
285db795 347void AliAnalysisTaskJetMatching::DoGeometricMatchingEtaPhi()
f1354962 348{
285db795 349 // do geometric matching based on eta phi distance between jet centers
90c6bce5 350 #ifdef DEBUGTASK
351 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
352 #endif
f1354962 353 fNoMatchedJets = 0; // reset the matched jet counter
354 Int_t iSource(fSourceJets->GetEntriesFast()), iTarget(fTargetJets->GetEntriesFast());
355 for(Int_t i(0); i < iSource; i++) {
356 AliEmcalJet* sourceJet(static_cast<AliEmcalJet*>(fSourceJets->At(i)));
4f6d96b5 357 if(!PassesCuts(sourceJet, 0)) continue;
f1354962 358 for(Int_t j(0); j < iTarget; j++) {
359 AliEmcalJet* targetJet(static_cast<AliEmcalJet*>(fTargetJets->At(j)));
4f6d96b5 360 if(!PassesCuts(targetJet, 1)) continue;
f6d1b1a7 361 if (fUseEmcalBaseJetCuts && !AcceptJet(targetJet, 1)) continue;
362 if((TMath::Abs(sourceJet->Eta() - targetJet->Eta()) < fMatchEta )) {
363 Double_t sourcePhi(sourceJet->Phi()), targetPhi(targetJet->Phi());
364 if(TMath::Abs(sourcePhi-targetPhi) > TMath::Abs(sourcePhi-targetPhi+TMath::TwoPi())) sourcePhi+=TMath::TwoPi();
365 if(TMath::Abs(sourcePhi-targetPhi) > TMath::Abs(sourcePhi-targetPhi-TMath::TwoPi())) sourcePhi-=TMath::TwoPi();
285db795 366 if(TMath::Abs(sourcePhi-targetPhi) < fMatchPhi) { // accept the jets as matching
367 Bool_t isBestMatch(kTRUE);
368 if(fGetBijection) { // match has been found, for bijection only store it if there's no better match
90c6bce5 369 #ifdef DEBUGTASK
370 printf(" > Entering first bijection test \n");
371 #endif
285db795 372 for(Int_t k(i); k < iSource; k++) {
373 AliEmcalJet* candidateSourceJet(static_cast<AliEmcalJet*>(fSourceJets->At(k)));
4f6d96b5 374 if(PassesCuts(candidateSourceJet, 0)) continue;
90c6bce5 375 #ifdef DEBUGTASK
376 printf("source distance %.2f \t candidate distance %.2f \n", GetR(sourceJet, targetJet),GetR(candidateSourceJet, targetJet));
377 #endif
285db795 378 if(GetR(sourceJet, targetJet) > GetR(candidateSourceJet, targetJet)) {
379 isBestMatch = kFALSE;
380 break;
381 }
382 }
90c6bce5 383 #ifdef DEBUGTASK
384 (isBestMatch) ? printf(" kept source \n ") : printf(" we can do better (rejected source) \n");
385 #endif
285db795 386 }
387 if(isBestMatch) {
388 fMatchedJetContainer[fNoMatchedJets][0] = sourceJet;
389 fMatchedJetContainer[fNoMatchedJets][1] = targetJet;
390 fNoMatchedJets++;
391 }
c2460f8a 392 if(fNoMatchedJets > 199) { // maximum to keep the cache at reasonable size
f6d1b1a7 393 AliError(Form("%s: Found too many matched jets (> 100). Adjust matching criteria !", GetName()));
394 return;
395 }
f1354962 396 }
397 }
398 }
399 }
400}
401//_____________________________________________________________________________
285db795 402void AliAnalysisTaskJetMatching::DoGeometricMatchingR()
f1354962 403{
404 // do geometric matching based on shortest path between jet centers
90c6bce5 405 #ifdef DEBUGTASK
406 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
407 #endif
f1354962 408 fNoMatchedJets = 0; // reset the matched jet counter
409 Int_t iSource(fSourceJets->GetEntriesFast()), iTarget(fTargetJets->GetEntriesFast());
410 for(Int_t i(0); i < iSource; i++) {
411 AliEmcalJet* sourceJet(static_cast<AliEmcalJet*>(fSourceJets->At(i)));
4f6d96b5 412 if(!PassesCuts(sourceJet, 0)) continue;
f1354962 413 for(Int_t j(0); j < iTarget; j++) {
414 AliEmcalJet* targetJet(static_cast<AliEmcalJet*>(fTargetJets->At(j)));
4f6d96b5 415 if(!PassesCuts(targetJet, 1)) continue;
f6d1b1a7 416 else if (fUseEmcalBaseJetCuts && !AcceptJet(targetJet, 1)) continue;
417 if(GetR(sourceJet, targetJet) <= fMatchR) {
285db795 418 Bool_t isBestMatch(kTRUE);
419 if(fGetBijection) { // match has been found, for bijection only store it if there's no better match
90c6bce5 420 #ifdef DEBUGTASK
421 printf(" > Entering first bijection test \n");
422 #endif
285db795 423 for(Int_t k(i); k < iSource; k++) {
424 AliEmcalJet* candidateSourceJet(static_cast<AliEmcalJet*>(fSourceJets->At(k)));
4f6d96b5 425 if(!PassesCuts(candidateSourceJet, 0)) continue;
90c6bce5 426 #ifdef DEBUGTASK
427 printf("source distance %.2f \t candidate distance %.2f \n", GetR(sourceJet, targetJet),GetR(candidateSourceJet, targetJet));
428 #endif
285db795 429 if(GetR(sourceJet, targetJet) > GetR(candidateSourceJet, targetJet)) {
430 isBestMatch = kFALSE;
431 break;
432 }
433 }
90c6bce5 434 #ifdef DEBUGTASK
435 (isBestMatch) ? printf(" kept source \n ") : printf(" we can do better (rejected source) \n");
436 #endif
285db795 437 }
438 if(isBestMatch) {
439 fMatchedJetContainer[fNoMatchedJets][0] = sourceJet;
440 fMatchedJetContainer[fNoMatchedJets][1] = targetJet;
441 fNoMatchedJets++;
442 }
f1354962 443 if(fNoMatchedJets > 99) {
444 AliError(Form("%s: Found too many matched jets (> 100). Adjust matching criteria !", GetName()));
445 return;
446 }
447 }
448 }
449 }
450}
451//_____________________________________________________________________________
08fae484 452void AliAnalysisTaskJetMatching::DoDiJetMatching()
453{
454 // match dijets. this is in a sense a 'special' mode of the task as both jet supplied jet arrays
455 // (target and source) are the same jet collection, matching will be performed on distribution in
456 // azimuth
26503cec 457 // no ouptut array is produced, dedicated dijet histo's are filled in this function instead
90c6bce5 458 #ifdef DEBUGTASK
459 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
460 #endif
26503cec 461 // get the leading jet in a given acceptance (TPC is default)
08fae484 462 Int_t leadingJetIndex(-1), subLeadingJetIndex(-1);
26503cec 463 // retrieve the leading jet, leadingJetIndex points to the leading jet
08fae484 464 AliEmcalJet* leadingJet(GetLeadingJet(fSourceJets, leadingJetIndex));
465 if(!leadingJet) return;
26503cec 466 // fill phi and eta of leading jet (should always be in selected acceptance)
10d7808d 467 fHistDiJetLeadingJet->Fill(leadingJet->Phi(), leadingJet->Eta(), leadingJet->Pt());
08fae484 468 Double_t sourcePhi(leadingJet->Phi()), targetPhi(-1);
26503cec 469 // get the sub-leading jet - faster when leading jet is also provided
08fae484 470 AliEmcalJet* subLeadingJet(GetSubLeadingJet(fSourceJets, leadingJetIndex, subLeadingJetIndex));
471 if(!subLeadingJet) return;
472 else { // check if the sub-leading jet is actually the away-side jet
473 targetPhi = subLeadingJet->Phi() + TMath::Pi();
26503cec 474 // rotate jets to common phase
475 if(TMath::Abs(sourcePhi) > TMath::Abs(sourcePhi+TMath::TwoPi())) sourcePhi+=TMath::TwoPi();
476 if(TMath::Abs(sourcePhi) > TMath::Abs(sourcePhi-TMath::TwoPi())) sourcePhi-=TMath::TwoPi();
477 if(TMath::Abs(targetPhi) > TMath::Abs(targetPhi+TMath::TwoPi())) targetPhi+=TMath::TwoPi();
478 if(TMath::Abs(targetPhi) > TMath::Abs(targetPhi-TMath::TwoPi())) targetPhi-=TMath::TwoPi();
479 if(TMath::Abs(sourcePhi-targetPhi) < fMatchPhi) {
10d7808d 480 fHistDiJet->Fill(subLeadingJet->Phi(), subLeadingJet->Eta(), leadingJet->Pt());
481 fHistDiJetDPhi->Fill(sourcePhi-targetPhi, leadingJet->Pt());
482 fHistDiJetDPt->Fill(leadingJet->Pt() - subLeadingJet->Pt(), leadingJet->Pt());
26503cec 483 }
08fae484 484 }
485}
486//_____________________________________________________________________________
285db795 487void AliAnalysisTaskJetMatching::DoConstituentMatching()
f1354962 488{
285db795 489 // match constituents within matched jets
90c6bce5 490 #ifdef DEBUGTASK
491 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
492 #endif
f1354962 493 if(!fTracks) {
494 AliFatal(Form("%s Fatal error! To do deep matching, supply jet constituents ! \n", GetName()));
495 return; // coverity ...
496 }
497 for(Int_t i(0); i < fNoMatchedJets; i++) {
0aaf31ab 498 AliEmcalJet* sourceJet = fMatchedJetContainer[i][0];
499 AliEmcalJet* targetJet = fMatchedJetContainer[i][1];
f1354962 500 if(sourceJet && targetJet) { // duplicate check: slot migth be NULL
0aaf31ab 501 Double_t targetPt(0);
502 Int_t iSJ(sourceJet->GetNumberOfTracks());
503 Int_t iTJ(targetJet->GetNumberOfTracks());
504 Int_t overlap(0), alreadyFound(0);
505 for(Int_t j(0); j < iSJ; j++) {
506 alreadyFound = 0;
507 Int_t idSource((Int_t)sourceJet->TrackAt(j));
508 for(Int_t k(0); k < iTJ; k++) { // compare all target tracks to the source track
509 if(idSource == targetJet->TrackAt(k) && alreadyFound == 0) {
510 overlap++;
511 alreadyFound++; // avoid possible duplicate matching
512 AliVParticle* vp(static_cast<AliVParticle*>(targetJet->TrackAt(k, fTracks)));
513 if(vp) targetPt += vp->Pt();
514 continue;
515 }
f1354962 516 }
517 }
285db795 518 if((float)overlap/(float)iSJ < fMinFracRecoveredConstituents || targetPt/sourceJet->Pt() < fMinFracRecoveredConstituentPt) {
90c6bce5 519 #ifdef DEBUGTASK
520 printf(" \n > Purging jet, recovered constituents ratio %i / %i = %.2f < or pt ratio %.2f / %.2f = %.2f < %.2f", overlap, iSJ, (float)overlap/(float)iSJ, targetPt, sourceJet->Pt(), targetPt/sourceJet->Pt(), fMinFracRecoveredConstituentPt);
521 #endif
f6d1b1a7 522 fMatchedJetContainer[i][0] = 0x0;
523 fMatchedJetContainer[i][1] = 0x0;
524 continue;
525 }
0aaf31ab 526 if(sourceJet->Pt() > 0) {
f6d1b1a7 527 Double_t sourceRho(0), targetRho(0);
528 if(fSourceRho) sourceRho = fSourceRho->GetLocalVal(sourceJet->Phi(), fSourceRadius)*sourceJet->Area();
529 if(fTargetRho) targetRho = fTargetRho->GetLocalVal(targetJet->Phi(), fTargetRadius)*targetJet->Area();
530 fProfFracPtMatched->Fill(sourceJet->Pt()-sourceRho, (targetPt-targetRho) / (sourceJet->Pt()-sourceRho));
531 fProfFracPtJets->Fill(sourceJet->Pt()-sourceRho, (targetJet->Pt()-targetRho) / (sourceJet->Pt()-sourceRho));
532 fProfFracNoMatched->Fill(sourceJet->Pt()-sourceRho, (double)overlap / (double)sourceJet->GetNumberOfTracks());
533 fProfFracNoJets->Fill(sourceJet->Pt()-sourceRho, (double)targetJet->GetNumberOfTracks() / (double)sourceJet->GetNumberOfTracks());
0aaf31ab 534 }
90c6bce5 535 #ifdef DEBUGTASK
536 if(fDebug > 0) {
537 printf("\n > Jet A: %i const\t", iSJ);
538 printf(" > Jet B %i const\t", iTJ);
539 printf(" -> OVERLAP: %i tracks <- \n", overlap);
f1354962 540 }
90c6bce5 541 #endif
542 }
f1354962 543 }
544}
545//_____________________________________________________________________________
285db795 546void AliAnalysisTaskJetMatching::GetBijection()
f1354962 547{
285db795 548 // bijection of source and matched jets, based on closest distance between jets
90c6bce5 549 #ifdef DEBUGTASK
550 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
551 #endif
f1354962 552 for(Int_t i(0); i < fNoMatchedJets; i++) {
553 for(Int_t j(i+1); j < fNoMatchedJets; j++) {
f6d1b1a7 554 if(fMatchedJetContainer[i][0] == fMatchedJetContainer[j][0]) {
555 // found source with two targets, now see which target is closer to the source
556 if(!(fMatchedJetContainer[i][0] && fMatchedJetContainer[i][1] && fMatchedJetContainer[j][0] && fMatchedJetContainer[j][1] )) continue;
285db795 557 Double_t rA(GetR(fMatchedJetContainer[i][0], fMatchedJetContainer[i][1])); // distance between connection A = i
558 Double_t rB(GetR(fMatchedJetContainer[j][0], fMatchedJetContainer[j][1])); // distance between connection B = j
559 if (rB > rA) { // jet two is far away, clear it from both target and source list
f6d1b1a7 560 fMatchedJetContainer[j][0] = 0x0;
561 fMatchedJetContainer[j][1] = 0x0;
285db795 562 } else { // jet one is far away, clear it from both target and source list
563 fMatchedJetContainer[i][0] = fMatchedJetContainer[j][0]; // switch to avoid breaking loop
564 fMatchedJetContainer[i][1] = fMatchedJetContainer[j][1];
565 fMatchedJetContainer[j][0] = 0x0; // clear duplicate jet from cache
f6d1b1a7 566 fMatchedJetContainer[j][1] = 0x0;
f1354962 567 }
90c6bce5 568 #ifdef DEBUGTASK
569 printf(" found duplicate jet, chose %.2f over %.2f \n" , (rB > rA) ? rA : rB, (rB > rA) ? rB : rA);
570 #endif
f1354962 571 }
572 }
573 }
574}
575//_____________________________________________________________________________
576void AliAnalysisTaskJetMatching::PostMatchedJets()
577{
578 // post matched jets
90c6bce5 579 #ifdef DEBUGTASK
580 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
581 #endif
459e1638 582 fMatchedJets->Delete(); // should be a NULL operation, but added just in case
f1354962 583 for(Int_t i(0), p(0); i < fNoMatchedJets; i++) {
584 if(fMatchedJetContainer[i][1]) { // duplicate jets will have NULL value here and are skipped
585 new((*fMatchedJets)[p]) AliEmcalJet(*fMatchedJetContainer[i][1]);
586 p++;
587 }
588 }
589}
590//_____________________________________________________________________________
591void AliAnalysisTaskJetMatching::FillAnalysisSummaryHistogram() const
592{
593 // fill the analysis summary histrogram, saves all relevant analysis settigns
90c6bce5 594 #ifdef DEBUGTASK
595 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
596 #endif
f1354962 597 fHistAnalysisSummary->GetXaxis()->SetBinLabel(1, "fUseScaledRho");
598 fHistAnalysisSummary->SetBinContent(1, (int)fUseScaledRho);
599 fHistAnalysisSummary->GetXaxis()->SetBinLabel(2, "fMatchingScheme");
600 fHistAnalysisSummary->SetBinContent(2, (int)fMatchingScheme);
601 fHistAnalysisSummary->GetXaxis()->SetBinLabel(3, "fSourceBKG");
602 fHistAnalysisSummary->SetBinContent(3, (int)fSourceBKG);
603 fHistAnalysisSummary->GetXaxis()->SetBinLabel(4, "fTargetBKG");
604 fHistAnalysisSummary->SetBinContent(4, (int)fTargetBKG);
285db795 605 fHistAnalysisSummary->GetXaxis()->SetBinLabel(5, "fMatchEta");
606 fHistAnalysisSummary->SetBinContent(5, fMatchEta);
607 fHistAnalysisSummary->GetXaxis()->SetBinLabel(6, "fMatchPhi");
608 fHistAnalysisSummary->SetBinContent(6, fMatchPhi);
609 fHistAnalysisSummary->GetXaxis()->SetBinLabel(7, "fMatchR");
610 fHistAnalysisSummary->SetBinContent(7, fMatchR);
611 fHistAnalysisSummary->GetXaxis()->SetBinLabel(8, "iter");
612 fHistAnalysisSummary->SetBinContent(8, 1.);
f1354962 613}
614//_____________________________________________________________________________
c2460f8a 615void AliAnalysisTaskJetMatching::FillMatchedJetHistograms()
f1354962 616{
617 // fill matched jet histograms
90c6bce5 618 #ifdef DEBUGTASK
619 printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
620 #endif
f1354962 621 for(Int_t i(0); i < fSourceJets->GetEntriesFast(); i++) {
622 AliEmcalJet* source = static_cast<AliEmcalJet*>(fSourceJets->At(i));
623 if(!source) continue;
c2460f8a 624 else if(fUseEmcalBaseJetCuts && !AcceptJet(source, 0)) continue;
f6d1b1a7 625 Double_t sourceRho(0), targetRho(0);
626 if(fSourceRho) sourceRho = fSourceRho->GetLocalVal(source->Phi(), fSourceRadius)*source->Area();
627 fHistSourceJetPt->Fill(source->Pt()-sourceRho);
f1354962 628 fHistNoConstSourceJet->Fill(source->GetNumberOfConstituents());
629 for(Int_t j(0); j < fTargetJets->GetEntriesFast(); j++) {
630 AliEmcalJet* target = static_cast<AliEmcalJet*>(fTargetJets->At(j));
631 if(target) {
c2460f8a 632 if(fUseEmcalBaseJetCuts && !AcceptJet(target, 1)) continue;
f6d1b1a7 633 if(fTargetRho) targetRho = fTargetRho->GetLocalVal(target->Phi(), fTargetRadius)*target->Area();
634 fProfQA->Fill(0.5, TMath::Abs((source->Pt()-sourceRho)-(target->Pt()-targetRho)));
635 fProfQA->Fill(1.5, TMath::Abs(source->Eta()-target->Eta()));
636 fProfQA->Fill(2.5, TMath::Abs(source->Phi()-target->Phi()));
f1354962 637 fHistUnsortedCorrelation->Fill(PhaseShift(source->Phi()-target->Phi(), 2));
638 if(j==0) {
f6d1b1a7 639 fHistTargetJetPt->Fill(target->Pt()-targetRho);
f1354962 640 fHistNoConstTargetJet->Fill(target->GetNumberOfConstituents());
641 }
642 }
643 }
644 }
645 for(Int_t i(0); i < fNoMatchedJets; i++) {
646 if(fMatchedJetContainer[i][0] && fMatchedJetContainer[i][1]) {
f6d1b1a7 647 Double_t sourceRho(0), targetRho(0);
648 if(fSourceRho) sourceRho = fSourceRho->GetLocalVal(fMatchedJetContainer[i][0]->Phi(), fSourceRadius)*fMatchedJetContainer[i][0]->Area();
968e3e0d 649 if(fTargetRho) targetRho = fTargetRho->GetLocalVal(fMatchedJetContainer[i][1]->Phi(), fTargetRadius)*fMatchedJetContainer[i][1]->Area();
f1354962 650 fHistMatchedCorrelation->Fill(PhaseShift(fMatchedJetContainer[i][0]->Phi()-fMatchedJetContainer[i][1]->Phi(), 2));
2e3ee924 651 fHistMatchedSourceJetPt->Fill(fMatchedJetContainer[i][0]->Pt()-sourceRho);
f6d1b1a7 652 fHistMatchedJetPt->Fill(fMatchedJetContainer[i][1]->Pt()-targetRho);
653 fHistNoConstMatchJet->Fill(fMatchedJetContainer[i][1]->Pt()-targetRho);
654 fProfQAMatched->Fill(0.5, TMath::Abs((fMatchedJetContainer[i][0]->Pt()-sourceRho)-(fMatchedJetContainer[i][1]->Pt()-targetRho)));
f1354962 655 fProfQAMatched->Fill(1.5, TMath::Abs(fMatchedJetContainer[i][0]->Eta()-fMatchedJetContainer[i][1]->Eta()));
656 fProfQAMatched->Fill(2.5, TMath::Abs(fMatchedJetContainer[i][0]->Phi()-fMatchedJetContainer[i][1]->Phi()));
4f6d96b5 657
c2460f8a 658 fHistSourceMatchedJetPt->Fill(fMatchedJetContainer[i][0]->Pt()-sourceRho, fMatchedJetContainer[i][1]->Pt()-targetRho);
2e3ee924 659 if(fDoDetectorResponse) {
660 fProfFracPtJets->Fill(fMatchedJetContainer[i][0]->Pt()-sourceRho, (fMatchedJetContainer[i][1]->Pt()-targetRho) / (fMatchedJetContainer[i][0]->Pt()-sourceRho));
661 fProfFracNoJets->Fill(fMatchedJetContainer[i][0]->Pt()-sourceRho, (double)fMatchedJetContainer[i][1]->GetNumberOfTracks() / (double)fMatchedJetContainer[i][0]->GetNumberOfTracks());
7ec64e26 662 fHistDetectorResponseProb->Fill((fMatchedJetContainer[i][1]->Pt()-fMatchedJetContainer[i][0]->Pt())/fMatchedJetContainer[i][0]->Pt(), fMatchedJetContainer[i][0]->Pt());
2e3ee924 663 }
f1354962 664 }
665 }
666}
667//_____________________________________________________________________________
08fae484 668AliEmcalJet* AliAnalysisTaskJetMatching::GetLeadingJet(TClonesArray* source, Int_t &leadingJetIndex, Double_t etaMin, Double_t etaMax)
669{
c52fd114 670 // return the leading jet within giiven acceptance
08fae484 671 Int_t iJets(source->GetEntriesFast());
672 Double_t pt(0);
673 AliEmcalJet* leadingJet(0x0);
674 for(Int_t i(0); i < iJets; i++) {
675 AliEmcalJet* jet = static_cast<AliEmcalJet*>(source->At(i));
676 if(jet->Eta() < etaMin || jet->Eta() > etaMax) continue;
677 if(jet->Pt() > pt) {
678 leadingJet = jet;
679 pt = leadingJet->Pt();
680 leadingJetIndex = i;
681 }
682 }
683 return leadingJet;
684}
685//_____________________________________________________________________________
c52fd114 686AliEmcalJet* AliAnalysisTaskJetMatching::GetSubLeadingJet(TClonesArray* source, Int_t leadingJetIndex, Int_t &subLeadingJetIndex)
08fae484 687{
c52fd114 688 // return the sub-leading jet within given acceptance
08fae484 689 // same as GetLeadingJet() but skips the leading jet (so returned jet is
690 // sub-leading by design)
c52fd114 691 // there is no eta requirement on the location of the sub-leading jet
08fae484 692 Int_t iJets(source->GetEntriesFast());
26503cec 693 // if the leading jet isn't given, retrieve it
c52fd114 694 if(leadingJetIndex < 0) GetLeadingJet(source, leadingJetIndex);
695 AliEmcalJet *leadingJet(0x0), *prevLeadingJet(static_cast<AliEmcalJet*>(source->At(leadingJetIndex)));
696 if(!prevLeadingJet) return 0x0;
697 Double_t pt(0), leadingPt(prevLeadingJet->Pt());
08fae484 698 for(Int_t i(0); i < iJets; i++) {
699 if(i == leadingJetIndex) continue;
700 AliEmcalJet* jet = static_cast<AliEmcalJet*>(source->At(i));
c52fd114 701 // check if jet is actually sub-leading
702 if(jet->Pt() > pt && jet->Pt() <= leadingPt) {
08fae484 703 leadingJet = jet;
704 pt = leadingJet->Pt();
705 subLeadingJetIndex = i;
706 }
707 }
708 return leadingJet;
709}
710//_____________________________________________________________________________
f1354962 711void AliAnalysisTaskJetMatching::PrintInfo() const
712{
713 // print some info
285db795 714 printf("\n > No. of source jets from %s \n \t %i \n ", fSourceJetsName.Data(), fSourceJets->GetEntriesFast());
f1354962 715 printf(" > No. of target jets from %s \n \t %i \n ", fTargetJetsName.Data(), fTargetJets->GetEntriesFast());
716 printf(" > No. of matched jets from %s \n \t %i \n ", fMatchedJetsName.Data(), fMatchedJets->GetEntriesFast());
f1354962 717}
718//_____________________________________________________________________________
719void AliAnalysisTaskJetMatching::Terminate(Option_t *)
720{
721 // terminate
722}
723//_____________________________________________________________________________