1 // **************************************
2 // Task used for estimating a charged to neutral correction
3 // sona.pochybova@cern.ch
4 // *******************************************
7 /**************************************************************************
8 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
10 * Author: The ALICE Off-line Project. *
11 * Contributors are mentioned in the code where appropriate. *
13 * Permission to use, copy, modify and distribute this software and its *
14 * documentation strictly for non-commercial purposes is hereby granted *
15 * without fee, provided that the above copyright notice appears in all *
16 * copies and that both the copyright notice and this permission notice *
17 * appear in the supporting documentation. The authors make no claims *
18 * about the suitability of this software for any purpose. It is *
19 * provided "as is" without express or implied warranty. *
20 **************************************************************************/
28 #include <TLorentzVector.h>
29 #include <TClonesArray.h>
30 #include <TRefArray.h>
34 #include "AliAnalysisTaskJetCorrections.h"
35 #include "AliAnalysisManager.h"
36 #include "AliAODEvent.h"
37 #include "AliAODVertex.h"
38 #include "AliAODHandler.h"
39 #include "AliAODTrack.h"
40 #include "AliAODJet.h"
41 //#include "AliMCEvent.h"
43 #include "AliAnalysisHelperJetTasks.h"
47 // corrections to jet energy by sona
51 ClassImp(AliAnalysisTaskJetCorrections)
53 AliAnalysisTaskJetCorrections::AliAnalysisTaskJetCorrections() : AliAnalysisTaskSE(),
95 for (Int_t i = 0; i < 3; i++)
100 fhECorrJet001[i] = 0;
104 fhdEvsErec001[i] = 0;
108 fhdPhidEta001[i] = 0;
109 fhdPhidEtaPt10[i] = 0;
110 fhdPhidEtaPt05[i] = 0;
111 fhdPhidEtaPt01[i] = 0;
112 fhdPhidEtaPt001[i] = 0;
116 AliAnalysisTaskJetCorrections::AliAnalysisTaskJetCorrections(const char * name):
117 AliAnalysisTaskSE(name),
124 fUseAODInput(kFALSE),
143 fhE2E1vsEsumGen(0x0),
144 fhE2E1vsEsumRec(0x0),
147 fhE2E1vsdPhiGen(0x0),
148 fhE2E1vsdPhiRec(0x0),
150 fhTrackBalance2(0x0),
151 fhTrackBalance3(0x0),
159 for (Int_t i = 0; i < 3; i++)
164 fhECorrJet001[i] = 0;
168 fhdEvsErec001[i] = 0;
172 fhdPhidEta001[i] = 0;
173 fhdPhidEtaPt10[i] = 0;
174 fhdPhidEtaPt05[i] = 0;
175 fhdPhidEtaPt01[i] = 0;
176 fhdPhidEtaPt001[i] = 0;
178 DefineOutput(1, TList::Class());
183 Bool_t AliAnalysisTaskJetCorrections::Notify()
186 // Implemented Notify() to read the cross sections
187 // and number of trials from pyxsec.root
191 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
192 Float_t xsection = 0;
195 Float_t fAvgTrials = 1;
197 TFile *curfile = tree->GetCurrentFile();
199 Error("Notify","No current file");
202 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
203 // construct a poor man average trials
204 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
205 if(ftrials>=nEntries)fAvgTrials = ftrials/nEntries; // CKB take this into account for normalisation
208 if(xsection>0)fXsection = xsection;
215 //___________________________________________________________________________________________________________________________________
216 void AliAnalysisTaskJetCorrections::UserCreateOutputObjects()
219 // Create the output container
221 // Printf("Analysing event %s :: # %5d\n", gSystem->pwd(), (Int_t) fEntry);
224 fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
226 Printf("%s:%d AODEvent not found in Input Manager %d",(char*)__FILE__,__LINE__,fUseAODInput);
231 // assume that the AOD is in the general output...
234 Printf("%s:%d AODEvent not found in the Output",(char*)__FILE__,__LINE__);
239 printf("AnalysisTaskJetSpectrum::UserCreateOutputObjects() \n");
243 fhEGen = new TH1F("EGen", "", 100, 0, 200);
247 fhERec = new TH1F("ERec", "", 100, 0, 200);
251 fhEGenRest = new TH1F("EGenRest", "", 100, 0, 200);
253 fList->Add(fhEGenRest);
255 fhERecRest = new TH1F("ERecRest", "", 100, 0, 200);
257 fList->Add(fhERecRest);
259 fhEsumGenRest = new TH1F("EsumGenRest", "", 100, 0, 200);
260 fhEsumGenRest->Sumw2();
261 fList->Add(fhEsumGenRest);
263 fhEsumRecRest = new TH1F("EsumRecRest", "", 100, 0, 200);
264 fhEsumRecRest->Sumw2();
265 fList->Add(fhEsumRecRest);
267 fhE2vsE1Gen = new TH2F("E2vsE1Gen", "", 100, 0, 200, 100, 0, 200);
268 fhE2vsE1Gen->Sumw2();
269 fList->Add(fhE2vsE1Gen);
271 fhE2vsE1Rec = new TH2F("E2vsE1Rec", "", 100, 0, 200, 100, 0, 200);
272 fhE2vsE1Rec->Sumw2();
273 fList->Add(fhE2vsE1Rec);
275 fhE2E1vsEsumGen = new TH2F("E2E1vsEsumGen", "", 100, 0, 200, 25, 0, 1);
276 fhE2E1vsEsumGen->Sumw2();
277 fList->Add(fhE2E1vsEsumGen);
279 fhE2E1vsEsumRec = new TH2F("E2E1vsEsumRec", "", 100, 0, 200, 25, 0, 1);
280 fhE2E1vsEsumRec->Sumw2();
281 fList->Add(fhE2E1vsEsumRec);
283 fhE2E1vsE1Gen = new TH2F("E2E1vsE1Gen", "", 100, 0, 200, 25, 0, 1);
284 fhE2E1vsE1Gen->Sumw2();
285 fList->Add(fhE2E1vsE1Gen);
287 fhE2E1vsE1Rec = new TH2F("E2E1vsE1Rec", "", 100, 0, 200, 25, 0, 1);
288 fhE2E1vsE1Rec->Sumw2();
289 fList->Add(fhE2E1vsE1Rec);
291 fhE2E1vsdPhiGen = new TH2F("E2E1vsdPhiGen", "", 64, -3.20, 3.20, 25, 0, 1);
292 fList->Add(fhE2E1vsdPhiGen);
294 fhE2E1vsdPhiRec = new TH2F("E2E1vsdPhiRec", "", 64, -3.20, 3.20, 25, 0, 1);
295 fList->Add(fhE2E1vsdPhiRec);
297 fhTrackBalance2 = new TH2F("TrackBalance2", "", 60, 0, 30, 60, 0, 30);
298 fhTrackBalance2->Sumw2();
299 fList->Add(fhTrackBalance2);
301 fhTrackBalance3 = new TH2F("TrackBalance3", "", 60, 0, 30, 60, 0, 30);
302 fhTrackBalance3->Sumw2();
303 fList->Add(fhTrackBalance3);
305 fhEt1Et22 = new TH2F("Et1Et22", "", 100, 0, 50, 100, 0, 50);
307 fList->Add(fhEt1Et22);
309 fhEt1Et23 = new TH2F("Et1Et23", "", 100, 0, 50, 100, 0, 50);
311 fList->Add(fhEt1Et23);
313 for(Int_t i = 0; i < 3; i++)
315 fhECorrJet10[i] = new TProfile(Form("ECorrJet10%d", i+1), "", 100, 0, 200, 0, 10);
316 fhECorrJet10[i]->SetXTitle("E_{rec} [GeV]");
317 fhECorrJet10[i]->SetYTitle("C=E_{gen}/E_{rec}");
318 fhECorrJet10[i]->Sumw2();
320 fhECorrJet05[i] = new TProfile(Form("ECorrJet05%d", i+1), "", 100, 0, 200, 0, 10);
321 fhECorrJet05[i]->SetXTitle("E_{rec} [GeV]");
322 fhECorrJet05[i]->SetYTitle("C=E_{gen}/E_{rec}");
323 fhECorrJet05[i]->Sumw2();
325 fhECorrJet01[i] = new TProfile(Form("ECorrJet01%d", i+1), "", 100, 0, 200, 0, 10);
326 fhECorrJet01[i]->SetXTitle("E_{rec} [GeV]");
327 fhECorrJet01[i]->SetYTitle("C=E_{gen}/E_{rec}");
328 fhECorrJet01[i]->Sumw2();
330 fhECorrJet001[i] = new TProfile(Form("ECorrJet001%d", i+1), "", 100, 0, 200, 0, 10);
331 fhECorrJet001[i]->SetXTitle("E_{rec} [GeV]");
332 fhECorrJet001[i]->SetYTitle("C=E_{gen}/E_{rec}");
333 fhECorrJet001[i]->Sumw2();
335 fhdEvsErec10[i] = new TProfile(Form("dEvsErec10_%d", i+1),"", 100, 0, 200, -1, 10);
336 fhdEvsErec10[i]->SetYTitle("|E_{rec}-E_{rec}|/E_{rec}");
337 fhdEvsErec10[i]->SetXTitle("E_{rec} [GeV]");
338 fhdEvsErec10[i]->Sumw2();
340 fhdEvsErec05[i] = new TProfile(Form("dEvsErec05_%d", i+1),"", 100, 0, 200, -1, 10);
341 fhdEvsErec05[i]->SetYTitle("|E_{rec}-E_{rec}|/E_{rec}");
342 fhdEvsErec05[i]->SetXTitle("E_{rec} [GeV]");
343 fhdEvsErec05[i]->Sumw2();
345 fhdEvsErec01[i] = new TProfile(Form("dEvsErec01_%d", i+1),"", 100, 0, 200, -1, 10);
346 fhdEvsErec01[i]->SetYTitle("|E_{rec}-E_{rec}|/E_{rec}");
347 fhdEvsErec01[i]->SetXTitle("E_{rec} [GeV]");
348 fhdEvsErec01[i]->Sumw2();
350 fhdEvsErec001[i] = new TProfile(Form("dEvsErec001_%d", i+1),"", 100, 0, 200, -1, 10);
351 fhdEvsErec001[i]->SetYTitle("|E_{rec}-E_{rec}|/E_{rec}");
352 fhdEvsErec001[i]->SetXTitle("E_{rec} [GeV]");
353 fhdEvsErec001[i]->Sumw2();
355 fhdPhidEta10[i] = new TH2F(Form("dPhidEta10_%d", i+1), "", 63, (-1)*TMath::Pi(), TMath::Pi(), 18, -0.9, 0.9);
356 fhdPhidEta10[i]->SetXTitle("#phi [rad]");
357 fhdPhidEta10[i]->SetYTitle("#eta");
358 fhdPhidEta10[i]->Sumw2();
360 fhdPhidEta05[i] = new TH2F(Form("dPhidEta05_%d", i+1), "", 63, (-1)*TMath::Pi(), TMath::Pi(), 18, -0.9, 0.9);
361 fhdPhidEta05[i]->SetXTitle("#phi [rad]");
362 fhdPhidEta05[i]->SetYTitle("#eta");
363 fhdPhidEta05[i]->Sumw2();
365 fhdPhidEta01[i] = new TH2F(Form("dPhidEta01_%d", i+1), "", 63, (-1)*TMath::Pi(), TMath::Pi(), 18, -0.9, 0.9);
366 fhdPhidEta01[i]->SetXTitle("#phi [rad]");
367 fhdPhidEta01[i]->SetYTitle("#eta");
368 fhdPhidEta01[i]->Sumw2();
370 fhdPhidEta001[i] = new TH2F(Form("dPhidEta001_%d", i+1), "", 63, (-1)*TMath::Pi(), TMath::Pi(), 18, -0.9, 0.9);
371 fhdPhidEta001[i]->SetXTitle("#phi [rad]");
372 fhdPhidEta001[i]->SetYTitle("#eta");
373 fhdPhidEta001[i]->Sumw2();
375 fhdPhidEtaPt10[i] = new TH2F(Form("dPhidEtaPt10_%d", i+1), "", 63, (-1)*TMath::Pi(), TMath::Pi(), 18, -0.9, 0.9);
376 fhdPhidEtaPt10[i]->SetXTitle("#phi [rad]");
377 fhdPhidEtaPt10[i]->SetYTitle("#eta");
378 fhdPhidEtaPt10[i]->Sumw2();
380 fhdPhidEtaPt05[i] = new TH2F(Form("dPhidEtaPt05_%d", i+1), "", 63, (-1)*TMath::Pi(), TMath::Pi(), 18, -0.9, 0.9);
381 fhdPhidEtaPt05[i]->SetXTitle("#phi [rad]");
382 fhdPhidEtaPt05[i]->SetYTitle("#eta");
383 fhdPhidEtaPt05[i]->Sumw2();
385 fhdPhidEtaPt01[i] = new TH2F(Form("dPhidEtaPt01_%d", i+1), "", 63, (-1)*TMath::Pi(), TMath::Pi(), 18, -0.9, 0.9);
386 fhdPhidEtaPt01[i]->SetXTitle("#phi [rad]");
387 fhdPhidEtaPt01[i]->SetYTitle("#eta");
388 fhdPhidEtaPt01[i]->Sumw2();
390 fhdPhidEtaPt001[i] = new TH2F(Form("dPhidEtaPt001_%d", i+1), "", 63, (-1)*TMath::Pi(), TMath::Pi(), 18, -0.9, 0.9);
391 fhdPhidEtaPt001[i]->SetXTitle("#phi [rad]");
392 fhdPhidEtaPt001[i]->SetYTitle("#eta");
393 fhdPhidEtaPt001[i]->Sumw2();
395 fList->Add(fhECorrJet10[i]);
396 fList->Add(fhECorrJet05[i]);
397 fList->Add(fhECorrJet01[i]);
398 fList->Add(fhECorrJet001[i]);
399 fList->Add(fhdEvsErec10[i]);
400 fList->Add(fhdEvsErec05[i]);
401 fList->Add(fhdEvsErec01[i]);
402 fList->Add(fhdEvsErec001[i]);
403 fList->Add(fhdPhidEta10[i]);
404 fList->Add(fhdPhidEta05[i]);
405 fList->Add(fhdPhidEta01[i]);
406 fList->Add(fhdPhidEta001[i]);
407 fList->Add(fhdPhidEtaPt10[i]);
408 fList->Add(fhdPhidEtaPt05[i]);
409 fList->Add(fhdPhidEtaPt01[i]);
410 fList->Add(fhdPhidEtaPt001[i]);
413 Printf("UserCreateOutputObjects finished\n");
416 //__________________________________________________________________________________________________________________________________________
417 void AliAnalysisTaskJetCorrections::Init()
419 printf("AliAnalysisJetCut::Init() \n");
422 //____________________________________________________________________________________________________________________________________________
423 void AliAnalysisTaskJetCorrections::UserExec(Option_t * )
425 // if (fDebug > 1) printf("Analysing event # %5d\n", (Int_t) fEntry);
428 //create an AOD handler
429 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
433 Printf("%s:%d no output aodHandler found Jet",(char*)__FILE__,__LINE__);
437 // AliMCEvent* mcEvent =MCEvent();
439 // Printf("%s:%d no mcEvent",(char*)__FILE__,__LINE__);
443 if (fDebug > 10)Printf("%s:%d",(char*)__FILE__,__LINE__);
446 AliAODVertex * pvtx = dynamic_cast<AliAODVertex*>(fAOD->GetPrimaryVertex());
449 AliAODJet genJets[kMaxJets];
452 AliAODJet recJets[kMaxJets];
455 //array of reconstructed jets from the AOD input
456 TClonesArray *aodRecJets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fBranchRec.Data()));
458 Printf("%s:%d no reconstructed Jet array with name %s in AOD",(char*)__FILE__,__LINE__,fBranchRec.Data());
462 // reconstructed jets
463 nRecJets = aodRecJets->GetEntries();
464 nRecJets = TMath::Min(nRecJets, kMaxJets);
466 for(int ir = 0;ir < nRecJets;++ir)
468 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodRecJets->At(ir));
473 // If we set a second branch for the input jets fetch this
474 TClonesArray * aodGenJets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fBranchGen.Data()));
478 printf("NO MC jets branch with name %s Found \n",fBranchGen.Data());
483 nGenJets = aodGenJets->GetEntries();
484 nGenJets = TMath::Min(nGenJets, kMaxJets);
486 for(Int_t ig =0 ; ig < nGenJets; ++ig)
488 AliAODJet * tmp = dynamic_cast<AliAODJet*>(aodGenJets->At(ig));
493 Double_t eRec[kMaxJets];
494 Double_t eGen[kMaxJets];
496 Double_t eRecRest[kMaxJets];
497 Double_t eGenRest[kMaxJets];
499 // AliAODJet jetRec[kMaxJets];
500 AliAODJet jetGen[kMaxJets];
502 Int_t idxRec[kMaxJets];
503 Int_t idxGen[kMaxJets];
505 // Double_t EsumRec = 0;
506 // Double_t EsumGen =0;
508 TLorentzVector vRec[kMaxJets];
509 TLorentzVector vGen[kMaxJets];
511 TLorentzVector vsumRec;
512 TLorentzVector vsumGen;
514 TVector3 pRec[kMaxJets];
515 TVector3 pGen[kMaxJets];
522 AliAODJet selJets[kMaxJets];
524 // loop for applying the separation cut
525 for(Int_t i = 0; i < nGenJets; i++)
529 selJets[nGenSel] = genJets[i];
536 for(Int_t j = 0; j < nGenJets; j++)
540 Double_t dRij = genJets[i].DeltaR(&genJets[j]);
542 if(dRij > 2*fR) tag++;
549 selJets[nGenSel] = genJets[i];
556 for (Int_t gj = 0; gj < nGenSel; gj++)
558 eGen[gj] = selJets[gj].E();
559 fhEGen->Fill(eGen[gj], fXsection);
562 TMath::Sort(nGenSel, eGen, idxGen);
563 for (Int_t ig = 0; ig < nGenSel; ig++)
564 jetGen[ig] = selJets[idxGen[ig]];
567 for (Int_t i = 0; i < nGenSel; ++i)
569 vGen[i].SetPxPyPzE(jetGen[i].Px(), jetGen[i].Py(), jetGen[i].Pz(), jetGen[i].E());
570 pGen[i].SetXYZ(vGen[i].Px(), vGen[i].Py(), vGen[i].Pz());
574 if(nGenSel > 1 && pGen[0].DeltaPhi(pGen[1]) > 2.8)
576 fhE2vsE1Gen->Fill(jetGen[0].E(), jetGen[1].E(), fXsection);
577 fhE2E1vsEsumGen->Fill(jetGen[0].E()+jetGen[1].E(), TMath::Abs(jetGen[0].E()-jetGen[1].E())/jetGen[0].E(), fXsection);
578 fhE2E1vsE1Gen->Fill(jetGen[0].E(), TMath::Abs(jetGen[0].E()-jetGen[1].E())/jetGen[0].E(), fXsection);
579 Double_t deltaPhi = (jetGen[0].Phi()-jetGen[1].Phi());
580 if(deltaPhi > TMath::Pi()) deltaPhi = deltaPhi - 2.*TMath::Pi();
581 if(deltaPhi < (-1.*TMath::Pi())) deltaPhi = deltaPhi + 2.*TMath::Pi();
582 fhE2E1vsdPhiGen->Fill(deltaPhi, TMath::Abs(jetGen[0].E()-jetGen[1].E())/jetGen[0].E(), fXsection);
585 Double_t fPxGen = vsumGen.Px();
586 Double_t fPyGen = vsumGen.Py();
587 Double_t fPzGen = vsumGen.Pz();
588 Double_t fEGen = vsumGen.E();
590 Double_t eSumGenRest = 0;
591 for (Int_t j = 0; j < nGenSel; j++)
593 vGen[j].Boost(-fPxGen/fEGen, -fPyGen/fEGen, -fPzGen/fEGen);
594 eGenRest[j] = vGen[j].E();
596 fhEGenRest->Fill(eGenRest[j], fXsection);
597 eSumGenRest += eGenRest[j];
601 fhEsumGenRest->Fill(eSumGenRest, fXsection);
603 //END VARIABLES FOR MC JETS ---------------
611 AliAODJet recSelJets[kMaxJets];
613 for(Int_t i = 0; i < nRecJets; i++)
617 recSelJets[nRecSel] = recJets[i];
624 for(Int_t j = 0; j < nRecJets; j++)
628 Double_t dRij = recJets[i].DeltaR(&recJets[j]);
630 if(dRij > 2*fR) tag1++;
635 if(tag1/counter1 == 1)
637 recSelJets[nRecSel] = recJets[i];
644 if(nRecSel == 0) return;
645 Printf("******NUMBER OF JETS AFTER DELTA R CUT : %d **********\n", nRecSel);
646 //sort rec/gen jets by energy in C.M.S
647 AliAODJet jetRecTmp[kMaxJets];
649 Double_t jetTrackPt[kTracks];
650 TLorentzVector jetTrackTmp[kTracks];
652 for (Int_t rj = 0; rj < nRecSel; rj++)
654 TRefArray * jetTracksAOD = dynamic_cast<TRefArray*>(recSelJets[rj].GetRefTracks());
655 if(!jetTracksAOD) continue;
656 if(jetTracksAOD->GetEntries() < 3) continue;
657 Int_t nJetTracks = 0;
658 for(Int_t j = 0; j < jetTracksAOD->GetEntries(); j++)
660 AliAODTrack * track = dynamic_cast<AliAODTrack*>(jetTracksAOD->At(j));
663 track->GetCovarianceXYZPxPyPz(cv);
664 if(cv[14] > 1000.) continue;
665 jetTrackPt[nTracks] = track->Pt();
666 jetTrackTmp[nTracks].SetPxPyPzE(track->Px(),track->Py(),track->Pz(),track->E());
670 if(nJetTracks < 4) continue;
671 jetRecTmp[nAccJets] = recSelJets[rj];
672 eRec[nAccJets] = recSelJets[rj].E();
673 fhERec->Fill(eRec[nAccJets], fXsection);
677 if(nAccJets == 0) return;
678 if(nTracks == 0) return;
680 Printf(" ************ Number of accepted jets : %d ************ \n", nAccJets);
682 AliAODJet jetRecAcc[kMaxJets];
683 TMath::Sort(nAccJets, eRec, idxRec);
684 for (Int_t rj = 0; rj < nAccJets; rj++)
685 jetRecAcc[rj] = jetRecTmp[idxRec[rj]];
687 //rest frame for reconstructed jets
688 for (Int_t i = 0; i < nAccJets; i++)
690 vRec[i].SetPxPyPzE(jetRecAcc[i].Px(), jetRecAcc[i].Py(), jetRecAcc[i].Pz(), jetRecAcc[i].E());
691 pRec[i].SetXYZ(vRec[i].Px(), vRec[i].Py(), vRec[i].Pz());
695 //check balance of two leading hadrons, deltaPhi > 2.
696 Int_t idxTrack[kTracks];
697 TMath::Sort(nTracks, jetTrackPt, idxTrack);
699 TLorentzVector jetTrack[kTracks];
700 for(Int_t iTr = 0; iTr < nTracks; iTr++)
701 jetTrack[iTr] = jetTrackTmp[idxTrack[iTr]];
704 while(jetTrack[0].DeltaPhi(jetTrack[n]) < 2.8)
709 for(Int_t iTr = 0; iTr < nTracks; iTr++)
711 if(TMath::Abs(jetTrack[0].DeltaPhi(jetTrack[iTr]) < 1.) && iTr != 0)
712 et1 += jetTrack[iTr].Et();
714 if(TMath::Abs(jetTrack[n].DeltaPhi(jetTrack[iTr]) < 1.) && iTr != n)
715 et2 += jetTrack[iTr].Et();
720 fhTrackBalance2->Fill(jetTrack[0].Et(), jetTrack[n].Et());
721 fhEt1Et22->Fill(et1, et2);
725 fhTrackBalance3->Fill(jetTrack[0].Et(), jetTrack[n].Et());
726 fhEt1Et23->Fill(et1, et2);
729 if(nAccJets > 1 && pRec[0].DeltaPhi(pRec[1]) > 2.8)
731 fhE2vsE1Rec->Fill(jetRecAcc[0].E(), jetRecAcc[1].E(), fXsection);
732 fhE2E1vsEsumRec->Fill(jetRecAcc[0].E()+jetRecAcc[1].E(), TMath::Abs(jetRecAcc[0].E()-jetRecAcc[1].E())/jetRecAcc[0].E(), fXsection);
733 fhE2E1vsE1Rec->Fill(jetRecAcc[0].E(), TMath::Abs(jetRecAcc[0].E()-jetRecAcc[1].E())/jetRecAcc[0].E(), fXsection);
734 Double_t deltaPhi = (jetRecAcc[0].Phi()-jetRecAcc[1].Phi());
735 if(deltaPhi > TMath::Pi()) deltaPhi = deltaPhi - 2.*TMath::Pi();
736 if(deltaPhi < (-1.*TMath::Pi())) deltaPhi = deltaPhi + 2.*TMath::Pi();
737 fhE2E1vsdPhiRec->Fill(-1*deltaPhi, TMath::Abs(jetRecAcc[0].E()-jetRecAcc[1].E())/jetRecAcc[0].E(), fXsection);
740 Double_t fPx = vsumRec.Px();
741 Double_t fPy = vsumRec.Py();
742 Double_t fPz = vsumRec.Pz();
743 Double_t fE = vsumRec.E();
745 Double_t eSumRecRest = 0;
746 for (Int_t j = 0; j < nAccJets; j++)
748 vRec[j].Boost(-fPx/fE, -fPy/fE, -fPz/fE);
749 eRecRest[j] = vRec[j].E();
751 fhERecRest->Fill(eRecRest[j], fXsection);
752 eSumRecRest += eRecRest[j];
755 fhEsumRecRest->Fill(eSumRecRest, fXsection);
758 Int_t iGenIndex[kMaxJets]; // Index of the generated jet for i-th rec -1 if none
759 Int_t iRecIndex[kMaxJets]; // Index of the rec jet for i-th gen -1 if none
761 for(int i = 0;i<kMaxJets;++i){
762 iGenIndex[i] = iRecIndex[i] = -1;
765 if(nAccJets == nGenSel)
767 AliAnalysisHelperJetTasks::GetClosestJets(jetGen,nGenSel,jetRecAcc,nAccJets,
768 iGenIndex,iRecIndex,0);
770 for(int ir = 0;ir < nAccJets;ir++){
771 Int_t ig = iGenIndex[ir];
772 if(ig>=0&&ig<nGenSel){
773 Double_t dPhi = TMath::Abs(jetRecAcc[ir].Phi()-jetGen[ig].Phi());
774 if(dPhi > TMath::Pi()) dPhi = dPhi - 2.*TMath::Pi();
775 if(dPhi < (-1.*TMath::Pi())) dPhi = dPhi + 2.*TMath::Pi();
776 Double_t sigma = TMath::Abs(jetGen[ig].E()-jetRecAcc[ir].E())/jetGen[ig].E();
777 dR = jetRecAcc[ir].DeltaR(&jetGen[ig]);
778 if(dR < 2*fR && dR >= fR)
784 fhdEvsErec10[0]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
785 fhECorrJet10[0]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
786 for(Int_t iTr = 0; iTr < nTracks; iTr++)
788 fhdPhidEtaPt10[0]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
789 fhdPhidEta10[0]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
795 fhdEvsErec10[1]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
796 fhECorrJet10[1]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
797 for(Int_t iTr = 0; iTr < nTracks; iTr++)
799 fhdPhidEtaPt10[1]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
800 fhdPhidEta10[1]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
806 fhdEvsErec10[2]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
807 fhECorrJet10[2]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
808 for(Int_t iTr = 0; iTr < nTracks; iTr++)
810 fhdPhidEtaPt10[2]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
811 fhdPhidEta10[2]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
817 if(dR < fR && dR >= 0.1)
823 fhdEvsErec05[0]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
824 fhECorrJet05[0]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
825 for(Int_t iTr = 0; iTr < nTracks; iTr++)
827 fhdPhidEtaPt05[0]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
828 fhdPhidEta05[0]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
834 fhdEvsErec05[1]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
835 fhECorrJet05[1]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
836 for(Int_t iTr = 0; iTr < nTracks; iTr++)
838 fhdPhidEtaPt05[1]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
839 fhdPhidEta05[1]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
845 fhdEvsErec05[2]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
846 fhECorrJet05[2]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
847 for(Int_t iTr = 0; iTr < nTracks; iTr++)
849 fhdPhidEtaPt05[2]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
850 fhdPhidEta05[2]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
856 if(dR < 0.1 && dR >= 0.01)
862 fhdEvsErec01[0]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
863 fhECorrJet01[0]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
864 for(Int_t iTr = 0; iTr < nTracks; iTr++)
866 fhdPhidEtaPt01[0]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
867 fhdPhidEta01[0]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
873 fhdEvsErec01[1]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
874 fhECorrJet01[1]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
875 for(Int_t iTr = 0; iTr < nTracks; iTr++)
877 fhdPhidEtaPt01[1]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
878 fhdPhidEta01[1]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
884 fhdEvsErec01[2]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
885 fhECorrJet01[2]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
886 for(Int_t iTr = 0; iTr < nTracks; iTr++)
888 fhdPhidEtaPt01[2]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
889 fhdPhidEta01[2]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
895 if(dR > 0.01) continue;
900 fhECorrJet001[0]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
901 fhdEvsErec001[0]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
902 for(Int_t iTr = 0; iTr < nTracks; iTr++)
904 fhdPhidEtaPt001[0]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
905 fhdPhidEta001[0]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
911 fhECorrJet001[1]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
912 fhdEvsErec001[1]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
913 for(Int_t iTr = 0; iTr < nTracks; iTr++)
915 fhdPhidEtaPt001[1]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
916 fhdPhidEta001[1]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
922 fhECorrJet001[2]->Fill(jetRecAcc[ir].E(), jetGen[ig].E()/jetRecAcc[ir].E(), fXsection);
923 fhdEvsErec001[2]->Fill(jetRecAcc[ir].E(), sigma, fXsection);
924 for(Int_t iTr = 0; iTr < nTracks; iTr++)
926 fhdPhidEtaPt001[2]->Fill(jetTrack[iTr].Phi(), jetTrack[iTr].Eta(), jetTrack[iTr].Pt());
927 fhdPhidEta001[2]->Fill(jetTrack[iTr].Phi(), jetTrack[ir].Eta());
934 // loop over reconstructed jets
937 Printf("%s:%d",(char*)__FILE__,__LINE__);
941 Printf("%s:%d Data Posted",(char*)__FILE__,__LINE__);
946 //__________________________________________________________________________________________________________________________________________________
947 void AliAnalysisTaskJetCorrections::Terminate(Option_t *)
949 printf("AnalysisJetCorrections::Terminate()");
953 //_______________________________________User defined functions_____________________________________________________________________________________