+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+
+//
+// task compares jets in two or three branches,
+// written for analysis of jet embedding in HI events
+//
+
+
#include "TChain.h"
#include "TTree.h"
#include "TMath.h"
#include "TH1F.h"
+#include "TH1I.h"
#include "TH2F.h"
#include "TH3F.h"
#include "THnSparse.h"
#include "AliAODEvent.h"
#include "AliAODJet.h"
+#include "AliAODHandler.h"
#include "AliAnalysisTaskJetResponseV2.h"
+using std::cout;
+using std::endl;
+
ClassImp(AliAnalysisTaskJetResponseV2)
AliAnalysisTaskJetResponseV2::AliAnalysisTaskJetResponseV2() :
AliAnalysisTaskSE(),
-fESD(0x0),
-fAOD(0x0),
-fBackgroundBranch(""),
-fIsPbPb(kTRUE),
-fOfflineTrgMask(AliVEvent::kAny),
-fMinContribVtx(1),
-fVtxZMin(-8.),
-fVtxZMax(8.),
-fEvtClassMin(0),
-fEvtClassMax(4),
-fCentMin(0.),
-fCentMax(100.),
-fNInputTracksMin(0),
-fNInputTracksMax(-1),
-fJetEtaMin(-.5),
-fJetEtaMax(.5),
-fJetPtMin(20.),
-fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
-fJetPtFractionMin(0.5),
-fNMatchJets(4),
-fMatchMaxDist(0.8),
-fKeepJets(kFALSE),
-fkNbranches(2),
-fkEvtClasses(12),
-fOutputList(0x0),
-fbEvent(kTRUE),
-fbJetsMismatch1(kTRUE),
-fbJetsMismatch2(kTRUE),
-fbJetsRp(kTRUE),
-fbJetsDeltaPt(kTRUE),
-fbJetsEta(kTRUE),
-fbJetsPhi(kTRUE),
-fbJetsArea(kTRUE),
-fbJetsBeforeCut1(kTRUE),
-fbJetsBeforeCut2(kTRUE),
-fHistEvtSelection(0x0),
-fHistJetSelection(0x0),
-fh2JetSelection(0x0),
-fhnEvent(0x0),
-fhnJetsMismatch1(0x0),
-fhnJetsMismatch2(0x0),
-fhnJetsRp(0x0),
-fhnJetsDeltaPt(0x0),
-fhnJetsEta(0x0),
-fhnJetsPhi(0x0),
-fhnJetsArea(0x0),
-fhnJetsBeforeCut1(0x0),
-fhnJetsBeforeCut2(0x0)
+ fESD(0x0),
+ fAOD(0x0),
+ fAODOut(0x0),
+ fAODExtension(0x0),
+ fNonStdFile(""),
+ fBackgroundBranch(""),
+ fIsPbPb(kTRUE),
+ fOfflineTrgMask(AliVEvent::kAny),
+ fMinContribVtx(1),
+ fVtxZMin(-8.),
+ fVtxZMax(8.),
+ fEvtClassMin(0),
+ fEvtClassMax(4),
+ fCentMin(0.),
+ fCentMax(100.),
+ fNInputTracksMin(0),
+ fNInputTracksMax(-1),
+ fJetEtaMin(-.5),
+ fJetEtaMax(.5),
+ fJetPtMin(20.),
+ fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
+ fJetPtFractionMin(0.5),
+ fNMatchJets(4),
+ fMatchMaxDist(0.8),
+ fKeepJets(kFALSE),
+ fkNbranches(2),
+ fkEvtClasses(12),
+ fOutputList(0x0),
+ fbEvent(kTRUE),
+ fbJetsMismatch1(kTRUE),
+ fbJetsMismatch2(kTRUE),
+ fbJetsRp(kTRUE),
+ fbJetsDeltaPt(kTRUE),
+ fbJetsEta(kTRUE),
+ fbJetsPhi(kTRUE),
+ fbJetsArea(kTRUE),
+ fbJets3Branches(kFALSE),
+ fbJetsBeforeCut1(kTRUE),
+ fbJetsBeforeCut2(kTRUE),
+ fHistEvtSelection(0x0),
+ fHistJetSelection(0x0),
+ fh2JetSelection(0x0),
+ fhnEvent(0x0),
+ fhnJetsMismatch1(0x0),
+ fhnJetsMismatch2(0x0),
+ fhnJetsRp(0x0),
+ fhnJetsDeltaPt(0x0),
+ fhnJetsEta(0x0),
+ fhnJetsPhi(0x0),
+ fhnJetsArea(0x0),
+ fhnJets3Branches(0x0),
+ fhnJetsBeforeCut1(0x0),
+ fhnJetsBeforeCut2(0x0)
{
- // default Constructor
+ // default Constructor
- fJetBranchName[0] = "";
- fJetBranchName[1] = "";
+ fJetBranchName[0] = "";
+ fJetBranchName[1] = "";
+ fJetBranchName[2] = "";
- fListJets[0] = new TList;
- fListJets[1] = new TList;
+ fListJets[0] = new TList;
+ fListJets[1] = new TList;
+ fListJets[2] = new TList;
}
AliAnalysisTaskJetResponseV2::AliAnalysisTaskJetResponseV2(const char *name) :
-AliAnalysisTaskSE(name),
-fESD(0x0),
-fAOD(0x0),
-fBackgroundBranch(""),
-fIsPbPb(kTRUE),
-fOfflineTrgMask(AliVEvent::kAny),
-fMinContribVtx(1),
-fVtxZMin(-8.),
-fVtxZMax(8.),
-fEvtClassMin(0),
-fEvtClassMax(4),
-fCentMin(0.),
-fCentMax(100.),
-fNInputTracksMin(0),
-fNInputTracksMax(-1),
-fJetEtaMin(-.5),
-fJetEtaMax(.5),
-fJetPtMin(20.),
-fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
-fJetPtFractionMin(0.5),
-fNMatchJets(4),
-fMatchMaxDist(0.8),
-fKeepJets(kFALSE),
-fkNbranches(2),
-fkEvtClasses(12),
-fOutputList(0x0),
-fbEvent(kTRUE),
-fbJetsMismatch1(kTRUE),
-fbJetsMismatch2(kTRUE),
-fbJetsRp(kTRUE),
-fbJetsDeltaPt(kTRUE),
-fbJetsEta(kTRUE),
-fbJetsPhi(kTRUE),
-fbJetsArea(kTRUE),
-fbJetsBeforeCut1(kTRUE),
-fbJetsBeforeCut2(kTRUE),
-fHistEvtSelection(0x0),
-fHistJetSelection(0x0),
-fh2JetSelection(0x0),
-fhnEvent(0x0),
-fhnJetsMismatch1(0x0),
-fhnJetsMismatch2(0x0),
-fhnJetsRp(0x0),
-fhnJetsDeltaPt(0x0),
-fhnJetsEta(0x0),
-fhnJetsPhi(0x0),
-fhnJetsArea(0x0),
-fhnJetsBeforeCut1(0x0),
-fhnJetsBeforeCut2(0x0)
+ AliAnalysisTaskSE(name),
+ fESD(0x0),
+ fAOD(0x0),
+ fAODOut(0x0),
+ fAODExtension(0x0),
+ fNonStdFile(""),
+ fBackgroundBranch(""),
+ fIsPbPb(kTRUE),
+ fOfflineTrgMask(AliVEvent::kAny),
+ fMinContribVtx(1),
+ fVtxZMin(-8.),
+ fVtxZMax(8.),
+ fEvtClassMin(0),
+ fEvtClassMax(4),
+ fCentMin(0.),
+ fCentMax(100.),
+ fNInputTracksMin(0),
+ fNInputTracksMax(-1),
+ fJetEtaMin(-.5),
+ fJetEtaMax(.5),
+ fJetPtMin(20.),
+ fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
+ fJetPtFractionMin(0.5),
+ fNMatchJets(4),
+ fMatchMaxDist(0.8),
+ fKeepJets(kFALSE),
+ fkNbranches(2),
+ fkEvtClasses(12),
+ fOutputList(0x0),
+ fbEvent(kTRUE),
+ fbJetsMismatch1(kTRUE),
+ fbJetsMismatch2(kTRUE),
+ fbJetsRp(kTRUE),
+ fbJetsDeltaPt(kTRUE),
+ fbJetsEta(kTRUE),
+ fbJetsPhi(kTRUE),
+ fbJetsArea(kTRUE),
+ fbJets3Branches(kFALSE),
+ fbJetsBeforeCut1(kTRUE),
+ fbJetsBeforeCut2(kTRUE),
+ fHistEvtSelection(0x0),
+ fHistJetSelection(0x0),
+ fh2JetSelection(0x0),
+ fhnEvent(0x0),
+ fhnJetsMismatch1(0x0),
+ fhnJetsMismatch2(0x0),
+ fhnJetsRp(0x0),
+ fhnJetsDeltaPt(0x0),
+ fhnJetsEta(0x0),
+ fhnJetsPhi(0x0),
+ fhnJetsArea(0x0),
+ fhnJets3Branches(0x0),
+ fhnJetsBeforeCut1(0x0),
+ fhnJetsBeforeCut2(0x0)
{
- // Constructor
+ // Constructor
- fJetBranchName[0] = "";
- fJetBranchName[1] = "";
+ fJetBranchName[0] = "";
+ fJetBranchName[1] = "";
+ fJetBranchName[2] = "";
- fListJets[0] = new TList;
- fListJets[1] = new TList;
+ fListJets[0] = new TList;
+ fListJets[1] = new TList;
+ fListJets[2] = new TList;
- DefineOutput(1, TList::Class());
+ DefineOutput(1, TList::Class());
}
AliAnalysisTaskJetResponseV2::~AliAnalysisTaskJetResponseV2()
{
- delete fListJets[0];
- delete fListJets[1];
+ delete fListJets[0];
+ delete fListJets[1];
+ delete fListJets[2];
}
-void AliAnalysisTaskJetResponseV2::SetBranchNames(const TString &branch1, const TString &branch2)
+void AliAnalysisTaskJetResponseV2::SetBranchNames(const TString &branch1, const TString &branch2, const TString &branch3)
{
- fJetBranchName[0] = branch1;
- fJetBranchName[1] = branch2;
+ fJetBranchName[0] = branch1;
+ fJetBranchName[1] = branch2;
+ fJetBranchName[2] = branch3;
+
+ if(strlen(fJetBranchName[2].Data()) ) {
+ fkNbranches = 3;
+ fbJets3Branches = kTRUE;
+ }
}
void AliAnalysisTaskJetResponseV2::Init()
{
- // check for jet branches
- if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
- AliError("Jet branch name not set.");
- }
+ // check for jet branches
+ if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
+ AliError("Jet branch name not set.");
+ }
}
void AliAnalysisTaskJetResponseV2::UserCreateOutputObjects()
{
- // Create histograms
- // Called once
- OpenFile(1);
- if(!fOutputList) fOutputList = new TList;
- fOutputList->SetOwner(kTRUE);
-
- Bool_t oldStatus = TH1::AddDirectoryStatus();
- TH1::AddDirectory(kFALSE);
-
-
- fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5);
- fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
- fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
- fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)");
- fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
- fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
- fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)");
-
- fHistJetSelection = new TH1I("fHistJetSelection", "jet selection", 8, -0.5, 7.5);
- fHistJetSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
- fHistJetSelection->GetXaxis()->SetBinLabel(2,"probes IN");
- fHistJetSelection->GetXaxis()->SetBinLabel(3,"no matching jet");
- fHistJetSelection->GetXaxis()->SetBinLabel(4,"not in list");
- fHistJetSelection->GetXaxis()->SetBinLabel(5,"fraction cut");
- fHistJetSelection->GetXaxis()->SetBinLabel(6,"acceptance cut");
- fHistJetSelection->GetXaxis()->SetBinLabel(7,"p_{T} cut");
- fHistJetSelection->GetXaxis()->SetBinLabel(8,"trigger exclude mask");
-
- fh2JetSelection = new TH2F("fh2JetSelection", "jet selection", 8, -0.5, 7.5,100,0.,200.);
- fh2JetSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
- fh2JetSelection->GetXaxis()->SetBinLabel(2,"probes IN");
- fh2JetSelection->GetXaxis()->SetBinLabel(3,"no matching jet");
- fh2JetSelection->GetXaxis()->SetBinLabel(4,"not in list");
- fh2JetSelection->GetXaxis()->SetBinLabel(5,"fraction cut");
- fh2JetSelection->GetXaxis()->SetBinLabel(6,"acceptance cut");
- fh2JetSelection->GetXaxis()->SetBinLabel(7,"p_{T} cut");
- fh2JetSelection->GetXaxis()->SetBinLabel(8,"trigger exclude mask");
-
-
- UInt_t entries = 0; // bit coded, see GetDimParams() below
- UInt_t opt = 0; // bit coded, default (0) or high resolution (1)
-
- if(fbEvent){
- entries = 1<<0 | 1<<1 | 1<<2 | 1<<26; // cent : nInpTrks : rp psi : pT hard bin
- opt = 1<<0 | 1<<1; // centrality and nInpTrks in high resolution
- fhnEvent = NewTHnSparseF("fhnEvent", entries, opt);
- }
+ // Create histograms
+ // Called once
+ OpenFile(1);
+ if(!fOutputList) fOutputList = new TList;
+ fOutputList->SetOwner(kTRUE);
+
+ Bool_t oldStatus = TH1::AddDirectoryStatus();
+ TH1::AddDirectory(kFALSE);
+
+
+ fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5);
+ fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
+ fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
+ fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)");
+ fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
+ fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
+ fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)");
+
+ fHistJetSelection = new TH1I("fHistJetSelection", "jet selection", 8, -0.5, 7.5);
+ fHistJetSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
+ fHistJetSelection->GetXaxis()->SetBinLabel(2,"probes IN");
+ fHistJetSelection->GetXaxis()->SetBinLabel(3,"no matching jet");
+ fHistJetSelection->GetXaxis()->SetBinLabel(4,"not in list");
+ fHistJetSelection->GetXaxis()->SetBinLabel(5,"fraction cut");
+ fHistJetSelection->GetXaxis()->SetBinLabel(6,"acceptance cut");
+ fHistJetSelection->GetXaxis()->SetBinLabel(7,"p_{T} cut");
+ fHistJetSelection->GetXaxis()->SetBinLabel(8,"trigger exclude mask");
+
+ fh2JetSelection = new TH2F("fh2JetSelection", "jet selection", 8, -0.5, 7.5,100,0.,200.);
+ fh2JetSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
+ fh2JetSelection->GetXaxis()->SetBinLabel(2,"probes IN");
+ fh2JetSelection->GetXaxis()->SetBinLabel(3,"no matching jet");
+ fh2JetSelection->GetXaxis()->SetBinLabel(4,"not in list");
+ fh2JetSelection->GetXaxis()->SetBinLabel(5,"fraction cut");
+ fh2JetSelection->GetXaxis()->SetBinLabel(6,"acceptance cut");
+ fh2JetSelection->GetXaxis()->SetBinLabel(7,"p_{T} cut");
+ fh2JetSelection->GetXaxis()->SetBinLabel(8,"trigger exclude mask");
+
+
+ UInt_t entries = 0; // bit coded, see GetDimParams() below
+ UInt_t opt = 0; // bit coded, default (0) or high resolution (1)
+
+ if(fbEvent){
+ entries = 1<<0 | 1<<1 | 1<<2 | 1<<26; // cent : nInpTrks : rp psi : pT hard bin
+ opt = 1<<0 | 1<<1; // centrality and nInpTrks in high resolution
+ fhnEvent = NewTHnSparseF("fhnEvent", entries, opt);
+ }
- if(fbJetsMismatch1){ // real mismatch (no related rec jet found)
- // cent : nInpTrks : rp bins : probe pt : probe eta : probe phi : probe area : pT hard bin
- entries = 1<<0 | 1<<1 | 1<<3 | 1<<6 | 1<<8 | 1<<10 | 1<<12 | 1<<26;
- opt = 1<<6 | 1<<8 | 1<<10;
- fhnJetsMismatch1 = NewTHnSparseF("fhnJetsMismatch1", entries, opt);
- }
-
- if(fbJetsMismatch2){ // acceptance + fraction cut
- // cent : nInpTrks : rp bins : jetPt(2x) : jetEta(2x) : deltaEta : deltaR : fraction : pT hard bin
- entries = 1<<0 | 1<<1 | 1<<3 | 1<<6 | 1<<7 | 1<<8 | 1<<9 | 1<<15 | 1<<17 | 1<<19 | 1<<26;
- opt = 1<<6 | 1<<7 | 1<<8 | 1<<9;
- fhnJetsMismatch2 = NewTHnSparseF("fhnJetsMismatch2", entries, opt);
-
- }
+ if(fbJetsMismatch1){ // real mismatch (no related rec jet found)
+ // cent : nInpTrks : rp bins : probe pt : probe eta : probe phi : probe area : pT hard bin
+ entries = 1<<0 | 1<<1 | 1<<3 | 1<<6 | 1<<8 | 1<<10 | 1<<12 | 1<<26;
+ opt = 1<<6 | 1<<8 | 1<<10;
+ fhnJetsMismatch1 = NewTHnSparseF("fhnJetsMismatch1", entries, opt);
+ }
+
+ if(fbJetsMismatch2){ // acceptance + fraction cut
+ // cent : nInpTrks : rp bins : jetPt(2x) : jetEta(2x) : deltaEta : deltaR : fraction : pT hard bin
+ entries = 1<<0 | 1<<1 | 1<<3 | 1<<6 | 1<<7 | 1<<8 | 1<<9 | 1<<15 | 1<<17 | 1<<19 | 1<<26;
+ opt = 1<<6 | 1<<7 | 1<<8 | 1<<9;
+ fhnJetsMismatch2 = NewTHnSparseF("fhnJetsMismatch2", entries, opt);
+
+ }
- if(fbJetsRp){
- // cent : nInpTrks : rp bins : rp wrt jet(2x) : probe pT : probe area: deltaPt : rp phi : rho : correction for RP : local rho
- /*
- entries = 1<<0 | 1<<1 | 1<<3 | 1<<4 | 1<<5 | 1<<6 | 1<<12 | 1<<14 | 1<<2 | 1<<22 | 1<<23 | 1<<24 | 1<<25;
- opt = 1<<4 | 1<<5;*/
- // cent : nInpTrks : rp bins : rp wrt jet(2x) : probe pT : deltaPt : rp phi : rho : correction for RP | pT hard bin
- entries = 1<<0 | 1<<1 | 1<<3 | 1<<4 | 1<<5 | 1<<6 | 1<<14 | 1<<2 | 1<<22 | 1<<23 | 1<<26;
- opt = 1<<4 | 1<<5;
- fhnJetsRp = NewTHnSparseF("fhnJetsRp", entries, opt);
- }
-
- // cent : nInpTrks : rp bins: deltaPt : jetPt(2x) : deltaArea : pT hard bin (hr delta pt)
- if(fbJetsDeltaPt){
- entries = 1<<0 | 1<<1 | 1<<3 | 1<<14 | 1<<6 | 1<<7 | 1<<18 | 1<<26;
- opt = 1<<1 | 1<<14 | 1<<6 | 1<<7 | 1<<18;
- fhnJetsDeltaPt = NewTHnSparseF("fhnJetsDeltaPt", entries, opt);
- }
-
- // cent : nInpTrks : rp bins : deltaPt : jetPt(2x) : deltaR : deltaEta : jetEta(2x) : pT hard bin (hr for eta)
- if(fbJetsEta){
- entries = 1<<0 | 1<<1 | 1<<3 | 1<<14 | 1<<6 | 1<<7 | 1<<17 | 1<<15 | 1<<8 | 1<<9 | 1<<26;
- opt = 1<<15 | 1<<8 | 1<<9;
- fhnJetsEta = NewTHnSparseF("fhnJetsEta", entries, opt);
- }
-
- // cent : nInpTrks : rp bins : jetPt(2x) : jetPhi(2x) : deltaPt : deltaPhi : pT hard bin
- if(fbJetsPhi){
- entries = 1<<0 | 1<<1 | 1<<3 | 1<<6 | 1<<7 | 1<<10 | 1<<11 | 1<<14 | 1<<16 | 1<<26;
- opt = 1<<10 | 1<<11;
- fhnJetsPhi = NewTHnSparseF("fhnJetsPhi", entries, opt);
- }
-
- // cent : nInpTrks : rp bins : deltaArea : jetArea(2x) : deltaR : fraction : distance next rec jet : pT next jet : deltaPt : jetPt(2x) : pT hard bin (hr for area)
- if(fbJetsArea){
- entries = 1<<0 | 1<<1 | 1<<3 | 1<<18 | 1<<12 | 1<<13 | 1<<17 | 1<<19 | 1<<20 | 1<<21 | 1<<14 | 1<<6 | 1<<7 | 1<<26;
- opt = 1<<18 | 1<<12 | 1<<13;
- fhnJetsArea = NewTHnSparseF("fhnJetsArea", entries, opt);
- }
-
-
- //before cut
-
- // cent : nInpTrks : rp bins : fraction : jetPt(2x) : jetEta(2x) : jetPhi(2x) (low resolution) (with fraction, eta, phi, pt cuts possible)
- if(fbJetsBeforeCut1){
- entries = 1<<0 | 1<<1 | 1<<3 | 1<<19 | 1<<6 | 1<<7 | 1<<8 | 1<<9 | 1<<10 | 1<<11;
- opt = 0;
- fhnJetsBeforeCut1 = NewTHnSparseF("fhnJetsBeforeCut1", entries, opt);
- }
-
- // cent : nInpTrks : rp bins : deltaPt : jetPt(2x) : deltaR : deltaEta : jetEta(2x) (low resolution)
- if(fbJetsBeforeCut2){
- entries = 1<<0 | 1<<1 | 1<<3 | 1<<14 | 1<<6 | 1<<7 | 1<<17 | 1<<15 | 1<<8 | 1<<9;
- opt = 0;
- fhnJetsBeforeCut2 = NewTHnSparseF("fhnJetsBeforeCut2", entries, opt);
- }
-
- fOutputList->Add(fHistEvtSelection);
- fOutputList->Add(fHistJetSelection);
- fOutputList->Add(fh2JetSelection);
- if(fbEvent) fOutputList->Add(fhnEvent);
- if(fbJetsMismatch1) fOutputList->Add(fhnJetsMismatch1);
- if(fbJetsMismatch2) fOutputList->Add(fhnJetsMismatch2);
- if(fbJetsRp) fOutputList->Add(fhnJetsRp);
- if(fbJetsDeltaPt) fOutputList->Add(fhnJetsDeltaPt);
- if(fbJetsEta) fOutputList->Add(fhnJetsEta);
- if(fbJetsPhi) fOutputList->Add(fhnJetsPhi);
- if(fbJetsArea) fOutputList->Add(fhnJetsArea);
- if(fbJetsBeforeCut1) fOutputList->Add(fhnJetsBeforeCut1);
- if(fbJetsBeforeCut2) fOutputList->Add(fhnJetsBeforeCut2);
-
- // =========== Switch on Sumw2 for all histos ===========
- for (Int_t i=0; i<fOutputList->GetEntries(); ++i) {
- TH1 *h1 = dynamic_cast<TH1*>(fOutputList->At(i));
- if (h1){
- h1->Sumw2();
- continue;
- }
- THnSparse *hn = dynamic_cast<THnSparse*>(fOutputList->At(i));
- if (hn){
- hn->Sumw2();
- }
- }
- TH1::AddDirectory(oldStatus);
-
- PostData(1, fOutputList);
+ if(fbJetsRp){
+ // cent : nInpTrks : rp bins : rp wrt jet(2x) : probe pT : probe area: deltaPt : rp phi : rho : correction for RP : local rho
+ /*
+ entries = 1<<0 | 1<<1 | 1<<3 | 1<<4 | 1<<5 | 1<<6 | 1<<12 | 1<<14 | 1<<2 | 1<<22 | 1<<23 | 1<<24 | 1<<25;
+ opt = 1<<4 | 1<<5;*/
+ // cent : nInpTrks : rp bins : rp wrt jet(2x) : probe pT : deltaPt : rp phi : rho : correction for RP | pT hard bin
+ entries = 1<<0 | 1<<1 | 1<<3 | 1<<4 | 1<<5 | 1<<6 | 1<<14 | 1<<2 | 1<<22 | 1<<23 | 1<<26;
+ opt = 1<<4 | 1<<5;
+ fhnJetsRp = NewTHnSparseF("fhnJetsRp", entries, opt);
+ }
+
+ // cent : nInpTrks : rp bins: deltaPt : jetPt(2x) : deltaArea : pT hard bin (hr delta pt)
+ if(fbJetsDeltaPt){
+ entries = 1<<0 | 1<<1 | 1<<3 | 1<<14 | 1<<6 | 1<<7 | 1<<18 | 1<<26;
+ opt = 1<<1 | 1<<14 | 1<<6 | 1<<7 | 1<<18;
+ fhnJetsDeltaPt = NewTHnSparseF("fhnJetsDeltaPt", entries, opt);
+ }
+
+ // cent : nInpTrks : rp bins : deltaPt : jetPt(2x) : deltaR : deltaEta : jetEta(2x) : pT hard bin (hr for eta)
+ if(fbJetsEta){
+ entries = 1<<0 | 1<<1 | 1<<3 | 1<<14 | 1<<6 | 1<<7 | 1<<17 | 1<<15 | 1<<8 | 1<<9 | 1<<26;
+ opt = 1<<15 | 1<<8 | 1<<9;
+ fhnJetsEta = NewTHnSparseF("fhnJetsEta", entries, opt);
+ }
+
+ // cent : nInpTrks : rp bins : jetPt(2x) : jetPhi(2x) : deltaPt : deltaPhi : pT hard bin
+ if(fbJetsPhi){
+ entries = 1<<0 | 1<<1 | 1<<3 | 1<<6 | 1<<7 | 1<<10 | 1<<11 | 1<<14 | 1<<16 | 1<<26;
+ opt = 1<<10 | 1<<11;
+ fhnJetsPhi = NewTHnSparseF("fhnJetsPhi", entries, opt);
+ }
+
+ // cent : nInpTrks : rp bins : deltaArea : jetArea(2x) : deltaR : fraction : distance next rec jet : pT next jet : deltaPt : jetPt(2x) : pT hard bin (hr for area)
+ if(fbJetsArea){
+ entries = 1<<0 | 1<<1 | 1<<3 | 1<<18 | 1<<12 | 1<<13 | 1<<17 | 1<<19 | 1<<20 | 1<<21 | 1<<14 | 1<<6 | 1<<7 | 1<<26;
+ opt = 1<<18 | 1<<12 | 1<<13;
+ fhnJetsArea = NewTHnSparseF("fhnJetsArea", entries, opt);
+ }
+
+ // cent : nInpTrks : jetPt(3x) : deltaPt : delta : jetArea(3x) : fraction(2x) : deltaR(1x) : pT hard bin
+ if(fbJets3Branches){
+ entries = 1<<0 | 1<<1 | 1<<6 | 1<<7 | 1<<27 | 1<<14 | 1<<28 | 1<<12 | 1<<13 | 1<<29 | 1<<19 | 1<<30 | 1<<17 | 1<<26;
+ opt = 1<<6 | 1<<7 | 1<<27 | 1<<14 | 1<<28;
+ fhnJets3Branches = NewTHnSparseF("fhnJets3Branches", entries, opt);
+ }
+
+
+
+
+ //before cut
+
+ // cent : nInpTrks : rp bins : fraction : jetPt(2x) : jetEta(2x) : jetPhi(2x) (low resolution) (with fraction, eta, phi, pt cuts possible)
+ if(fbJetsBeforeCut1){
+ entries = 1<<0 | 1<<1 | 1<<3 | 1<<19 | 1<<6 | 1<<7 | 1<<8 | 1<<9 | 1<<10 | 1<<11;
+ opt = 0;
+ fhnJetsBeforeCut1 = NewTHnSparseF("fhnJetsBeforeCut1", entries, opt);
+ }
+
+ // cent : nInpTrks : rp bins : deltaPt : jetPt(2x) : deltaR : deltaEta : jetEta(2x) (low resolution)
+ if(fbJetsBeforeCut2){
+ entries = 1<<0 | 1<<1 | 1<<3 | 1<<14 | 1<<6 | 1<<7 | 1<<17 | 1<<15 | 1<<8 | 1<<9;
+ opt = 0;
+ fhnJetsBeforeCut2 = NewTHnSparseF("fhnJetsBeforeCut2", entries, opt);
+ }
+
+ fOutputList->Add(fHistEvtSelection);
+ fOutputList->Add(fHistJetSelection);
+ fOutputList->Add(fh2JetSelection);
+ if(fbEvent) fOutputList->Add(fhnEvent);
+ if(fbJetsMismatch1) fOutputList->Add(fhnJetsMismatch1);
+ if(fbJetsMismatch2) fOutputList->Add(fhnJetsMismatch2);
+ if(fbJetsRp) fOutputList->Add(fhnJetsRp);
+ if(fbJetsDeltaPt) fOutputList->Add(fhnJetsDeltaPt);
+ if(fbJetsEta) fOutputList->Add(fhnJetsEta);
+ if(fbJetsPhi) fOutputList->Add(fhnJetsPhi);
+ if(fbJetsArea) fOutputList->Add(fhnJetsArea);
+ if(fbJets3Branches) fOutputList->Add(fhnJets3Branches);
+ if(fbJetsBeforeCut1) fOutputList->Add(fhnJetsBeforeCut1);
+ if(fbJetsBeforeCut2) fOutputList->Add(fhnJetsBeforeCut2);
+
+ // =========== Switch on Sumw2 for all histos ===========
+ for (Int_t i=0; i<fOutputList->GetEntries(); ++i) {
+ TH1 *h1 = dynamic_cast<TH1*>(fOutputList->At(i));
+ if (h1){
+ h1->Sumw2();
+ continue;
+ }
+ THnSparse *hn = dynamic_cast<THnSparse*>(fOutputList->At(i));
+ if (hn){
+ hn->Sumw2();
+ }
+ }
+ TH1::AddDirectory(oldStatus);
+
+ PostData(1, fOutputList);
}
void AliAnalysisTaskJetResponseV2::UserExec(Option_t *)
{
- // load events, apply event cuts, then compare leading jets
-
- if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
- AliError("Jet branch name not set.");
- return;
- }
-
- fESD=dynamic_cast<AliESDEvent*>(InputEvent());
- if (!fESD) {
- AliError("ESD not available");
- fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
- } else {
- fAOD = dynamic_cast<AliAODEvent*>(AODEvent());
- }
- if (!fAOD) {
- AliError("AOD not available");
- return;
- }
-
- // -- event selection --
- fHistEvtSelection->Fill(1); // number of events before event selection
-
- // physics selection
- AliInputEventHandler* inputHandler = (AliInputEventHandler*)
- ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
- if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
- if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
- fHistEvtSelection->Fill(2);
- PostData(1, fOutputList);
- return;
- }
-
- // vertex selection
- AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
- Int_t nTracksPrim = primVtx->GetNContributors();
- if ((nTracksPrim < fMinContribVtx) ||
- (primVtx->GetZ() < fVtxZMin) ||
- (primVtx->GetZ() > fVtxZMax) ){
- if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
- fHistEvtSelection->Fill(3);
- PostData(1, fOutputList);
- return;
- }
-
- // event class selection (from jet helper task)
- Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
- if(fDebug) Printf("Event class %d", eventClass);
- if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
- fHistEvtSelection->Fill(4);
- PostData(1, fOutputList);
- return;
- }
-
- // centrality selection
- AliCentrality *cent = 0x0;
- Float_t centValue = 0.;
- if(fESD) cent = fESD->GetCentrality();
- if(cent) centValue = cent->GetCentralityPercentile("V0M");
- if(fDebug) printf("centrality: %f\n", centValue);
- if (centValue < fCentMin || centValue > fCentMax){
- fHistEvtSelection->Fill(4);
- PostData(1, fOutputList);
- return;
- }
-
-
- // multiplicity due to input tracks
- Int_t nInputTracks = GetNInputTracks();
-
- if (nInputTracks < fNInputTracksMin || (fNInputTracksMax > -1 && nInputTracks > fNInputTracksMax)){
- fHistEvtSelection->Fill(5);
- PostData(1, fOutputList);
- return;
- }
+ // load events, apply event cuts, then compare leading jets
+
+ if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
+ AliError("Jet branch name not set.");
+ return;
+ }
+
+ fESD=dynamic_cast<AliESDEvent*>(InputEvent());
+ if (!fESD) {
+ AliError("ESD not available");
+
+ fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
+ // assume that the AOD is in the general output...
+ fAODOut = AODEvent();
+ } else {
+ fAOD = dynamic_cast<AliAODEvent*>(AODEvent());
+ }
+ if (!fAOD) {
+ AliError("AOD not available");
+ return;
+ }
+
+ if(fNonStdFile.Length()!=0){
+ // case that we have an AOD extension we need can fetch the jets from the extended output
+ AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
+ fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
+ if(!fAODExtension){
+ if(fDebug>1)Printf("AODExtension found for %s",fNonStdFile.Data());
+ }
+ }
+
+ // -- event selection --
+ fHistEvtSelection->Fill(1); // number of events before event selection
+
+ // physics selection
+ AliInputEventHandler* inputHandler = (AliInputEventHandler*)
+ ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
+ if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
+ if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
+ fHistEvtSelection->Fill(2);
+ PostData(1, fOutputList);
+ return;
+ }
+
+ // vertex selection
+ AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
+ Int_t nTracksPrim = primVtx->GetNContributors();
+ if ((nTracksPrim < fMinContribVtx) ||
+ (primVtx->GetZ() < fVtxZMin) ||
+ (primVtx->GetZ() > fVtxZMax) ){
+ if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
+ fHistEvtSelection->Fill(3);
+ PostData(1, fOutputList);
+ return;
+ }
+
+ // event class selection (from jet helper task)
+ Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
+ if(fDebug) Printf("Event class %d", eventClass);
+ if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
+ fHistEvtSelection->Fill(4);
+ PostData(1, fOutputList);
+ return;
+ }
+
+ // centrality selection
+ AliCentrality *cent = 0x0;
+ Float_t centValue = 0.;
+ if(fESD) cent = fESD->GetCentrality();
+ if(cent) centValue = cent->GetCentralityPercentile("V0M");
+ if(fDebug) printf("centrality: %f\n", centValue);
+ if (centValue < fCentMin || centValue > fCentMax){
+ fHistEvtSelection->Fill(4);
+ PostData(1, fOutputList);
+ return;
+ }
+
+
+ // multiplicity due to input tracks
+ Int_t nInputTracks = GetNInputTracks();
+
+ if (nInputTracks < fNInputTracksMin || (fNInputTracksMax > -1 && nInputTracks > fNInputTracksMax)){
+ fHistEvtSelection->Fill(5);
+ PostData(1, fOutputList);
+ return;
+ }
- fHistEvtSelection->Fill(0); // accepted events
- // -- end event selection --
-
- // pt hard
- Double_t pthard = AliAnalysisTaskFastEmbedding::GetPtHard();
- Int_t pthardbin = GetPtHardBin(pthard);
-
- // reaction plane
- Double_t rp = AliAnalysisHelperJetTasks::ReactionPlane(kFALSE);
- if(fbEvent){
- Double_t eventEntries[3] = { (Double_t)centValue, (Double_t)nInputTracks, rp };
- fhnEvent->Fill(eventEntries);
- }
-
-
- // get background
- AliAODJetEventBackground* externalBackground = 0;
- if(!externalBackground&&fBackgroundBranch.Length()){
- externalBackground = (AliAODJetEventBackground*)(fAOD->FindListObject(fBackgroundBranch.Data()));
- //if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
- }
- Float_t rho = 0;
- if(externalBackground)rho = externalBackground->GetBackground(0);
-
-
- // fetch jets
- TClonesArray *aodJets[2];
- aodJets[0] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[0].Data())); // in general: embedded jet
- aodJets[1] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[1].Data())); // in general: embedded jet + UE
-
-
-
- for (Int_t iJetType = 0; iJetType < 2; iJetType++) {
- fListJets[iJetType]->Clear();
- if (!aodJets[iJetType]) continue;
-
- if(fDebug) Printf("%s: %d jets",fJetBranchName[iJetType].Data(),aodJets[iJetType]->GetEntriesFast());
+ fHistEvtSelection->Fill(0); // accepted events
+ // -- end event selection --
+
+ // pt hard
+ Double_t pthard = AliAnalysisTaskFastEmbedding::GetPtHard();
+ Int_t pthardbin = GetPtHardBin(pthard);
+
+ // reaction plane
+ Double_t rp = AliAnalysisHelperJetTasks::ReactionPlane(kFALSE);
+ if(fbEvent){
+ Double_t eventEntries[3] = { (Double_t)centValue, (Double_t)nInputTracks, rp };
+ fhnEvent->Fill(eventEntries);
+ }
+
+
+ // get background
+ AliAODJetEventBackground* externalBackground = 0;
+ if(!externalBackground&&fBackgroundBranch.Length()){
+ externalBackground = (AliAODJetEventBackground*)(fAOD->FindListObject(fBackgroundBranch.Data()));
+ if(!externalBackground && fAODOut) externalBackground = (AliAODJetEventBackground*)(fAODOut->FindListObject(fBackgroundBranch.Data()));
+ if(!externalBackground && fAODExtension) externalBackground = (AliAODJetEventBackground*)(fAODExtension->GetAOD()->FindListObject(fBackgroundBranch.Data()));
+ //if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
+ }
+ Float_t rho = 0;
+ if(externalBackground) rho = externalBackground->GetBackground(0);
+
+
+ // fetch jets
+ TClonesArray *aodJets[3];
+ aodJets[0] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[0].Data())); // in general: embedded jet
+ if(!aodJets[0] && fAODOut) aodJets[0] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[0].Data())); // in general: embedded jet
+ if(!aodJets[0] && fAODExtension) aodJets[0] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[0].Data()));
+ aodJets[1] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[1].Data())); // in general: embedded jet + UE version1
+ if(!aodJets[1] && fAODOut) aodJets[1] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[1].Data())); // in general: embedded jet
+ if(!aodJets[1] && fAODExtension) aodJets[1] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[1].Data()));
+ if( strlen(fJetBranchName[2].Data()) ) {
+ aodJets[2] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[2].Data())); // in general: embedded jet + UE version2
+ if(!aodJets[2] && fAODOut) aodJets[2] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[2].Data())); // in general: embedded jet
+ if(!aodJets[2] && fAODExtension) aodJets[2] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[2].Data()));
+ if(aodJets[2]) fkNbranches=3;
+ if(fDebug>10) printf("3rd branch: %s\n",fJetBranchName[2].Data());
+ }
+
+ if(fDebug>10) printf("fkNbranches %d\n",fkNbranches);
+ for (Int_t iJetType = 0; iJetType < fkNbranches; iJetType++) {
+ fListJets[iJetType]->Clear();
+ if (!aodJets[iJetType]) {
+ if(fDebug) Printf("%s: no jet branch\n",fJetBranchName[iJetType].Data());
+ continue;
+ }
+ if(fDebug) Printf("%s: %d jets",fJetBranchName[iJetType].Data(),aodJets[iJetType]->GetEntriesFast());
- for (Int_t iJet = 0; iJet < aodJets[iJetType]->GetEntriesFast(); iJet++) {
- AliAODJet *jet = dynamic_cast<AliAODJet*>((*aodJets[iJetType])[iJet]);
- if (jet) fListJets[iJetType]->Add(jet);
- }
- fListJets[iJetType]->Sort();
- }
+ for (Int_t iJet = 0; iJet < aodJets[iJetType]->GetEntriesFast(); iJet++) {
+ AliAODJet *jet = dynamic_cast<AliAODJet*>((*aodJets[iJetType])[iJet]);
+ if (jet) fListJets[iJetType]->Add(jet);
+ }
+ fListJets[iJetType]->Sort();
+ }
- // jet matching
- static TArrayI aMatchIndex(fListJets[0]->GetEntries());
- static TArrayF aPtFraction(fListJets[0]->GetEntries());
- if(aMatchIndex.GetSize()<fListJets[0]->GetEntries()) aMatchIndex.Set(fListJets[0]->GetEntries());
- if(aPtFraction.GetSize()<fListJets[0]->GetEntries()) aPtFraction.Set(fListJets[0]->GetEntries());
-
- // stores matched jets in 'aMatchIndex' and fraction of pT in 'aPtFraction'
- AliAnalysisHelperJetTasks::GetJetMatching(fListJets[0], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[0]->GetEntries()),
- fListJets[1], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[1]->GetEntries()),
- aMatchIndex, aPtFraction, fDebug, fMatchMaxDist, fIsPbPb?1:2);
+ // jet matching
+ static TArrayI aMatchIndex(fListJets[0]->GetEntries());
+ static TArrayF aPtFraction(fListJets[0]->GetEntries());
+ if(aMatchIndex.GetSize()<fListJets[0]->GetEntries()) aMatchIndex.Set(fListJets[0]->GetEntries());
+ if(aPtFraction.GetSize()<fListJets[0]->GetEntries()) aPtFraction.Set(fListJets[0]->GetEntries());
+
+ // stores matched jets in 'aMatchIndex' and fraction of pT in 'aPtFraction'
+ AliAnalysisHelperJetTasks::GetJetMatching(fListJets[0], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[0]->GetEntries()),
+ fListJets[1], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[1]->GetEntries()),
+ aMatchIndex, aPtFraction, fDebug, fMatchMaxDist, fIsPbPb?1:2);
+ static TArrayI aMatchIndexv2(fListJets[0]->GetEntries());
+ static TArrayF aPtFractionv2(fListJets[0]->GetEntries());
+ if(aMatchIndexv2.GetSize()<fListJets[0]->GetEntries()) aMatchIndexv2.Set(fListJets[0]->GetEntries());
+ if(aPtFractionv2.GetSize()<fListJets[0]->GetEntries()) aPtFractionv2.Set(fListJets[0]->GetEntries());
+ if( strlen(fJetBranchName[2].Data()) ) {
+ AliAnalysisHelperJetTasks::GetJetMatching(fListJets[0], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[0]->GetEntries()),
+ fListJets[2], TMath::Min((Int_t)fNMatchJets,(Int_t)fListJets[2]->GetEntries()),
+ aMatchIndexv2, aPtFractionv2, fDebug, fMatchMaxDist, fIsPbPb?1:2);
+ }
- // loop over matched jets
- Int_t ir = -1; // index of matched reconstruced jet
- Float_t fraction = -1.;
- AliAODJet *jet[2] = { 0x0, 0x0 };
- Float_t jetEta[2] = { -990., -990. };
- Float_t jetPhi[2] = { -990., -990. };
- Float_t jetPt[2] = { -990., -990. };
- Float_t jetArea[2] = { -990., -990. };
- Float_t rpJet[2] = { -990., -990. };
- Int_t rpBin = -990;
+ // loop over matched jets
+ Int_t ir = -1; // index of matched reconstruced jet
+ Int_t ir2 = -1; // index of matched reconstruced jet
+ Float_t fraction = -1.;
+ Float_t fraction2 = -1.;
+ AliAODJet *jet[3] = { 0x0, 0x0, 0x0 };
+ Float_t jetEta[3] = { -990., -990., -990. };
+ Float_t jetPhi[3] = { -990., -990., -990. };
+ Float_t jetPt[3] = { -990., -990., -990. };
+ Float_t jetArea[3] = { -990., -990., -990. };
+ Float_t rpJet[3] = { -990., -990., -990. };
+ Int_t rpBin = -990;
- for(Int_t ig=0; ig<fListJets[0]->GetEntries(); ++ig){
- ir = aMatchIndex[ig];
-
- //fetch jets
- jet[0] = (AliAODJet*)(fListJets[0]->At(ig));
- if(ir>=0) jet[1] = (AliAODJet*)(fListJets[1]->At(ir));
- else jet[1] = 0x0;
-
- for(Int_t i=0; i<fkNbranches; ++i){
- if(!jet[i]){
- jetEta[i] = -990;
- jetPhi[i] = -990.;
- jetPt[i] = -990.;
- jetArea[i] = -990.;
- rpJet[i] = -990.;
- if(i==1) rpBin = -990;
- }
- else {
- jetEta[i] = jet[i]->Eta();
- jetPhi[i] = jet[i]->Phi();
- jetPt[i] = GetPt(jet[i], i);
- jetArea[i] = jet[i]->EffectiveAreaCharged();
- rpJet[i] = TVector2::Phi_mpi_pi(rp-jetPhi[i]);
- if(i==1) rpBin = AliAnalysisHelperJetTasks::GetPhiBin(TVector2::Phi_mpi_pi(rp-jetPhi[i]), 3);
- }
+ for(Int_t ig=0; ig<fListJets[0]->GetEntries(); ++ig){
+ ir = aMatchIndex[ig];
+ if(aMatchIndexv2.GetSize()>ig) ir2 = aMatchIndexv2[ig];
+ else ir2 =0;
+
+ //fetch jets
+ jet[0] = (AliAODJet*)(fListJets[0]->At(ig));
+ if(ir>=0) jet[1] = (AliAODJet*)(fListJets[1]->At(ir));
+ else jet[1] = 0x0;
+ if(ir2>=0) jet[2] = (AliAODJet*)(fListJets[2]->At(ir2));
+ else jet[2] = 0x0;
+
+ for(Int_t i=0; i<fkNbranches; ++i){
+ if(!jet[i]){
+ jetEta[i] = -990;
+ jetPhi[i] = -990.;
+ jetPt[i] = -990.;
+ jetArea[i] = -990.;
+ rpJet[i] = -990.;
+ if(i==1) rpBin = -990;
+ }
+ else {
+ jetEta[i] = jet[i]->Eta();
+ jetPhi[i] = jet[i]->Phi();
+ jetPt[i] = GetPt(jet[i], i);
+ jetArea[i] = jet[i]->EffectiveAreaCharged();
+ rpJet[i] = TVector2::Phi_mpi_pi(rp-jetPhi[i]);
+ if(i==1) rpBin = AliAnalysisHelperJetTasks::GetPhiBin(TVector2::Phi_mpi_pi(rp-jetPhi[i]), 3);
}
- fraction = aPtFraction[ig];
-
- // jet statistics
- fHistJetSelection->Fill(1); // all probe jets
- if(jet[0]) fh2JetSelection->Fill(1.,jetPt[0]); // all probe jets
-
- if(ir<0){
- fHistJetSelection->Fill(2);
- if(jet[0]) fh2JetSelection->Fill(2.,jetPt[0]);
+ }
+ fraction = aPtFraction[ig];
+ if(aPtFractionv2.GetSize()>ig) fraction2 = aPtFractionv2[ig];
+ else fraction2 = 0.;
+
+ // jet statistics
+ fHistJetSelection->Fill(1); // all probe jets
+ if(jet[0]) fh2JetSelection->Fill(1.,jetPt[0]); // all probe jets
+
+ if(ir<0){
+ fHistJetSelection->Fill(2);
+ if(jet[0]) fh2JetSelection->Fill(2.,jetPt[0]);
- if(fbJetsMismatch1){
- if(!jet[0]) continue;
- Double_t jetEntriesMismatch1[7] = {
- (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
- (Double_t)jetPt[0], (Double_t)jetEta[0], (Double_t)jetPhi[0], (Double_t)jetArea[0]
- };
- fhnJetsMismatch1->Fill(jetEntriesMismatch1);
- }
-
- continue;
+ if(fbJetsMismatch1){
+ if(!jet[0]) continue;
+ Double_t jetEntriesMismatch1[7] = {
+ (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
+ (Double_t)jetPt[0], (Double_t)jetEta[0], (Double_t)jetPhi[0], (Double_t)jetArea[0]
+ };
+ fhnJetsMismatch1->Fill(jetEntriesMismatch1);
}
+
+ continue;
+ }
- if(!jet[0] || !jet[1]){
- fHistJetSelection->Fill(3);
- if(jet[0]) fh2JetSelection->Fill(3.,jetPt[0]);
- continue;
- }
+ if(!jet[0] || !jet[1]){
+ fHistJetSelection->Fill(3);
+ if(jet[0]) fh2JetSelection->Fill(3.,jetPt[0]);
+ continue;
+ }
- // look for distance to next rec jet
- Float_t distNextJet = -0.01; // no neighbor
- Float_t ptNextJet = -1.;
- for(Int_t r=0; r<fListJets[1]->GetEntries(); ++r){
- if(r==ir) continue;
- Float_t tmpDeltaR = jet[1]->DeltaR((AliAODJet*)fListJets[1]->At(r));
- if(distNextJet<0. || distNextJet>tmpDeltaR){
- distNextJet = tmpDeltaR;
- if(fKeepJets) ptNextJet = ((AliAODJet*)fListJets[1]->At(r))->GetPtSubtracted(0);
- else ptNextJet = ((AliAODJet*)fListJets[1]->At(r))->Pt();
- }
+ // look for distance to next rec jet
+ Float_t distNextJet = -0.01; // no neighbor
+ Float_t ptNextJet = -1.;
+ for(Int_t r=0; r<fListJets[1]->GetEntries(); ++r){
+ if(r==ir) continue;
+ Float_t tmpDeltaR = jet[1]->DeltaR((AliAODJet*)fListJets[1]->At(r));
+ if(distNextJet<0. || distNextJet>tmpDeltaR){
+ distNextJet = tmpDeltaR;
+ if(fKeepJets) ptNextJet = ((AliAODJet*)fListJets[1]->At(r))->GetPtSubtracted(0);
+ else ptNextJet = ((AliAODJet*)fListJets[1]->At(r))->Pt();
}
+ }
- //Float_t localRho = jetArea[1]>0. ? (jetPt[1]+rho*jetArea[1] - jetPt[0]) / jetArea[1] : 0.;
- //Float_t relRho = rho>0. ? localRho / rho : 0.;
+ //Float_t localRho = jetArea[1]>0. ? (jetPt[1]+rho*jetArea[1] - jetPt[0]) / jetArea[1] : 0.;
+ //Float_t relRho = rho>0. ? localRho / rho : 0.;
- // calculate parameters of associated jets
- /* from Leticia, kT clusterizer
- Float_t par0[4] = { 0.00409, 0.01229, 0.05, 0.26 };
- Float_t par1[4] = { -2.97035e-03, -2.03182e-03, -1.25702e-03, -9.95107e-04 };
- Float_t par2[4] = { 1.02865e-01, 1.49039e-01, 1.53910e-01, 1.51109e-01 };
- */
- // own, from embedded tracks
- Float_t par0[4] = { 0.02841, 0.05039, 0.09092, 0.24089 };
- Float_t par1[4] = { -4.26725e-04, -1.15273e-03, -1.56827e-03, -3.08003e-03 };
- Float_t par2[4] = { 4.95415e-02, 9.79538e-02, 1.32814e-01, 1.71743e-01 };
+ // calculate parameters of associated jets
+ /* from Leticia, kT clusterizer
+ Float_t par0[4] = { 0.00409, 0.01229, 0.05, 0.26 };
+ Float_t par1[4] = { -2.97035e-03, -2.03182e-03, -1.25702e-03, -9.95107e-04 };
+ Float_t par2[4] = { 1.02865e-01, 1.49039e-01, 1.53910e-01, 1.51109e-01 };
+ */
+ // own, from embedded tracks
+ Float_t par0[4] = { 0.02841, 0.05039, 0.09092, 0.24089 };
+ Float_t par1[4] = { -4.26725e-04, -1.15273e-03, -1.56827e-03, -3.08003e-03 };
+ Float_t par2[4] = { 4.95415e-02, 9.79538e-02, 1.32814e-01, 1.71743e-01 };
- Float_t rpCorr = 0.;
+ Float_t rpCorr = 0.;
- if(eventClass>0&&eventClass<4){
- rpCorr = par0[eventClass-1] + par1[eventClass-1] * 2*TMath::Cos(rpJet[1]) + par2[eventClass-1] * 2*TMath::Cos(2*rpJet[1]);
- rpCorr *= rho * jetArea[1];
- }
-
- Float_t deltaPt = jetPt[1]-jetPt[0];
- Float_t deltaEta = jetEta[1]-jetEta[0];
- Float_t deltaPhi = TVector2::Phi_mpi_pi(jetPhi[1]-jetPhi[0]);
- Float_t deltaR = TMath::Sqrt(deltaEta*deltaEta+deltaPhi*deltaPhi);
- Float_t deltaArea = jetArea[1]-jetArea[0];
+ if(eventClass>0&&eventClass<4){
+ rpCorr = par0[eventClass-1] + par1[eventClass-1] * 2*TMath::Cos(rpJet[1]) + par2[eventClass-1] * 2*TMath::Cos(2*rpJet[1]);
+ rpCorr *= rho * jetArea[1];
+ }
+
+ Float_t delta = jetPt[2]-jetPt[1];
+ Float_t deltaPt = jetPt[1]-jetPt[0];
+ Float_t deltaEta = jetEta[1]-jetEta[0];
+ Float_t deltaPhi = TVector2::Phi_mpi_pi(jetPhi[1]-jetPhi[0]);
+ Float_t deltaR = TMath::Sqrt(deltaEta*deltaEta+deltaPhi*deltaPhi);
+ Float_t deltaArea = jetArea[1]-jetArea[0];
- // fill thnsparse before acceptance cut
- if(fbJetsBeforeCut1){
- Double_t jetBeforeCutEntries1[10] = {
- (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
- (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetEta[0], (Double_t)jetEta[1],
- (Double_t)jetPhi[0], (Double_t)jetPhi[1], (Double_t)fraction
- };
- fhnJetsBeforeCut1->Fill(jetBeforeCutEntries1);
- }
+ // fill thnsparse before acceptance cut
+ if(fbJetsBeforeCut1){
+ Double_t jetBeforeCutEntries1[10] = {
+ (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
+ (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetEta[0], (Double_t)jetEta[1],
+ (Double_t)jetPhi[0], (Double_t)jetPhi[1], (Double_t)fraction
+ };
+ fhnJetsBeforeCut1->Fill(jetBeforeCutEntries1);
+ }
- if(fbJetsBeforeCut2){
- Double_t jetBeforeCutEntries2[10] = {
- (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
- (Double_t)jetPt[0], (Double_t)jetPt[1],
- (Double_t)jetEta[0], (Double_t)jetEta[1],
- (Double_t)deltaPt, (Double_t)deltaEta, (Double_t)deltaR
- };
- fhnJetsBeforeCut2->Fill(jetBeforeCutEntries2);
- }
+ if(fbJetsBeforeCut2){
+ Double_t jetBeforeCutEntries2[10] = {
+ (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
+ (Double_t)jetPt[0], (Double_t)jetPt[1],
+ (Double_t)jetEta[0], (Double_t)jetEta[1],
+ (Double_t)deltaPt, (Double_t)deltaEta, (Double_t)deltaR
+ };
+ fhnJetsBeforeCut2->Fill(jetBeforeCutEntries2);
+ }
- // jet selection
- Bool_t jetAccepted = kTRUE;
- // minimum fraction required
- if(fraction<fJetPtFractionMin){
- fHistJetSelection->Fill(4);
- fh2JetSelection->Fill(4.,jetPt[0]);
- jetAccepted = kFALSE;
- }
+ // jet selection
+ Bool_t jetAccepted = kTRUE;
+ // minimum fraction required
+ if(fraction<fJetPtFractionMin){
+ fHistJetSelection->Fill(4);
+ fh2JetSelection->Fill(4.,jetPt[0]);
+ jetAccepted = kFALSE;
+ }
- if(jetAccepted){
- // jet acceptance + minimum pT check
- if(jetEta[0]>fJetEtaMax || jetEta[0]<fJetEtaMin ||
- jetEta[1]>fJetEtaMax || jetEta[1]<fJetEtaMin){
+ if(jetAccepted){
+ // jet acceptance + minimum pT check
+ if(jetEta[0]>fJetEtaMax || jetEta[0]<fJetEtaMin ||
+ jetEta[1]>fJetEtaMax || jetEta[1]<fJetEtaMin){
- if(fDebug){
- Printf("Jet not in eta acceptance.");
- Printf("[0]: jet %d eta %.2f", ig, jetEta[0]);
- Printf("[1]: jet %d eta %.2f", ir, jetEta[1]);
- }
- fHistJetSelection->Fill(5);
- fh2JetSelection->Fill(5.,jetPt[0]);
- jetAccepted = kFALSE;
- }
+ if(fDebug){
+ Printf("Jet not in eta acceptance.");
+ Printf("[0]: jet %d eta %.2f", ig, jetEta[0]);
+ Printf("[1]: jet %d eta %.2f", ir, jetEta[1]);
+ }
+ fHistJetSelection->Fill(5);
+ fh2JetSelection->Fill(5.,jetPt[0]);
+ jetAccepted = kFALSE;
}
- if(jetAccepted){
- if(jetPt[1] < fJetPtMin){
- if(fDebug) Printf("Jet %d (pT %.1f GeV/c) has less than required pT.", ir, jetPt[1]);
- fHistJetSelection->Fill(6);
- fh2JetSelection->Fill(6.,jetPt[0]);
- jetAccepted = kFALSE;
- }
+ }
+ if(jetAccepted){
+ if(jetPt[0] < fJetPtMin){
+ if(fDebug) Printf("Jet %d (pT %.1f GeV/c) has less than required pT.", ir, jetPt[0]);
+ fHistJetSelection->Fill(6);
+ fh2JetSelection->Fill(6.,jetPt[0]);
+ jetAccepted = kFALSE;
}
+ }
- if(jetAccepted){
- if(jet[1]->Trigger()&fJetTriggerExcludeMask){
- fHistJetSelection->Fill(7);
- fh2JetSelection->Fill(7.,jetPt[0]);
- jetAccepted = kFALSE;
- }
-
+ if(jetAccepted){
+ if(jet[1]->Trigger()&fJetTriggerExcludeMask){
+ fHistJetSelection->Fill(7);
+ fh2JetSelection->Fill(7.,jetPt[0]);
+ jetAccepted = kFALSE;
}
+
+ }
- if(!jetAccepted){
- if(fbJetsMismatch2){
- Double_t jetEntriesMismatch2[10] = {
- (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
- (Double_t)jetPt[0], (Double_t)jetPt[1],
- (Double_t)jetEta[0], (Double_t)jetEta[1],
- (Double_t)deltaEta, (Double_t)deltaR,
- (Double_t)fraction
- };
- fhnJetsMismatch2->Fill(jetEntriesMismatch2);
- }
- continue;
+ if(!jetAccepted){
+ if(fbJetsMismatch2){
+ Double_t jetEntriesMismatch2[10] = {
+ (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
+ (Double_t)jetPt[0], (Double_t)jetPt[1],
+ (Double_t)jetEta[0], (Double_t)jetEta[1],
+ (Double_t)deltaEta, (Double_t)deltaR,
+ (Double_t)fraction
+ };
+ fhnJetsMismatch2->Fill(jetEntriesMismatch2);
}
+ continue;
+ }
- // all accepted jets
- fHistJetSelection->Fill(0);
- fh2JetSelection->Fill(0.,jetPt[0]);
+ // all accepted jets
+ fHistJetSelection->Fill(0);
+ fh2JetSelection->Fill(0.,jetPt[0]);
- // fill thnsparse
- if(fbJetsRp){
- Double_t jetEntriesRp[11] = {
- (Double_t)centValue, (Double_t)nInputTracks, (Double_t) rp,
- (Double_t)rpBin, (Double_t)rpJet[0], (Double_t)rpJet[1],
- (Double_t)jetPt[0], (Double_t)deltaPt, (Double_t)rho, (Double_t)rpCorr,
- (Double_t)pthardbin
- };
- fhnJetsRp->Fill(jetEntriesRp);
- }
+ // fill thnsparse
+ if(fbJetsRp){
+ Double_t jetEntriesRp[11] = {
+ (Double_t)centValue, (Double_t)nInputTracks, (Double_t) rp,
+ (Double_t)rpBin, (Double_t)rpJet[0], (Double_t)rpJet[1],
+ (Double_t)jetPt[0], (Double_t)deltaPt, (Double_t)rho, (Double_t)rpCorr,
+ (Double_t)pthardbin
+ };
+ fhnJetsRp->Fill(jetEntriesRp);
+ }
- if(fbJetsDeltaPt){
- Double_t jetEntriesDeltaPt[8] = {
- (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
- (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)deltaPt, (Double_t)deltaArea,
- (Double_t)pthardbin
- };
- fhnJetsDeltaPt->Fill(jetEntriesDeltaPt);
- }
+ if(fbJetsDeltaPt){
+ Double_t jetEntriesDeltaPt[8] = {
+ (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
+ (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)deltaPt, (Double_t)deltaArea,
+ (Double_t)pthardbin
+ };
+ fhnJetsDeltaPt->Fill(jetEntriesDeltaPt);
+ }
- if(fbJetsEta){
- Double_t jetEntriesEta[11] = {
- (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
- (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetEta[0], (Double_t)jetEta[1],
- (Double_t)deltaPt, (Double_t)deltaEta, (Double_t)deltaR,
- (Double_t)pthardbin
- };
- fhnJetsEta->Fill(jetEntriesEta);
- }
+ if(fbJetsEta){
+ Double_t jetEntriesEta[11] = {
+ (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
+ (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetEta[0], (Double_t)jetEta[1],
+ (Double_t)deltaPt, (Double_t)deltaEta, (Double_t)deltaR,
+ (Double_t)pthardbin
+ };
+ fhnJetsEta->Fill(jetEntriesEta);
+ }
- if(fbJetsPhi){
- Double_t jetEntriesPhi[10] = {
- (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
- (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetPhi[0], (Double_t)jetPhi[1],
- (Double_t)deltaPt, (Double_t)deltaPhi,
- (Double_t)pthardbin
- };
- fhnJetsPhi->Fill(jetEntriesPhi);
- }
+ if(fbJetsPhi){
+ Double_t jetEntriesPhi[10] = {
+ (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
+ (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetPhi[0], (Double_t)jetPhi[1],
+ (Double_t)deltaPt, (Double_t)deltaPhi,
+ (Double_t)pthardbin
+ };
+ fhnJetsPhi->Fill(jetEntriesPhi);
+ }
- if(fbJetsArea){
- Double_t jetEntriesArea[14] = {
- (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
- (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetArea[0], (Double_t)jetArea[1],
- (Double_t)deltaPt, (Double_t)deltaR, (Double_t)deltaArea,
- (Double_t)fraction, (Double_t)distNextJet, (Double_t)ptNextJet,
- (Double_t)pthardbin
- };
- fhnJetsArea->Fill(jetEntriesArea);
- }
+ if(fbJetsArea){
+ Double_t jetEntriesArea[14] = {
+ (Double_t)centValue, (Double_t)nInputTracks, (Double_t)rpBin,
+ (Double_t)jetPt[0], (Double_t)jetPt[1], (Double_t)jetArea[0], (Double_t)jetArea[1],
+ (Double_t)deltaPt, (Double_t)deltaR, (Double_t)deltaArea,
+ (Double_t)fraction, (Double_t)distNextJet, (Double_t)ptNextJet,
+ (Double_t)pthardbin
+ };
+ fhnJetsArea->Fill(jetEntriesArea);
+ }
+
+ if(fbJets3Branches){
+ Double_t jetEntries3Branches[14] = {
+ (Double_t)centValue, (Double_t)nInputTracks,
+ (Double_t)jetPt[0], (Double_t)jetPt[1],
+ (Double_t)jetArea[0], (Double_t)jetArea[1],
+ (Double_t)deltaPt, (Double_t)fraction, (Double_t)pthardbin,
+ (Double_t)jetPt[2],(Double_t)delta,(Double_t)jetArea[2], (Double_t)fraction2, (Double_t)deltaR
+ };
+ fhnJets3Branches->Fill(jetEntries3Branches);
+ }
- }
+ }
- PostData(1, fOutputList);
+ PostData(1, fOutputList);
}
void AliAnalysisTaskJetResponseV2::Terminate(const Option_t *)
{
- // Draw result to the screen
- // Called once at the end of the query
+ // Draw result to the screen
+ // Called once at the end of the query
- if (!GetOutputData(1))
- return;
+ if (!GetOutputData(1))
+ return;
}
Int_t AliAnalysisTaskJetResponseV2::GetNInputTracks()
{
-
- Int_t nInputTracks = 0;
-
- TString jbname(fJetBranchName[1]);
- //needs complete event, use jets without background subtraction
- for(Int_t i=1; i<=3; ++i){
- if(jbname.Contains(Form("B%d",i))) jbname.ReplaceAll(Form("B%d",i),"B0");
- }
- // use only HI event
- if(jbname.Contains("AODextraonly")) jbname.ReplaceAll("AODextraonly","AOD");
- if(jbname.Contains("AODextra")) jbname.ReplaceAll("AODextra","AOD");
-
- if(fDebug) Printf("Multiplicity from jet branch %s", jbname.Data());
- TClonesArray *tmpAODjets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(jbname.Data()));
- if(!tmpAODjets){
- Printf("Jet branch %s not found", jbname.Data());
- Printf("AliAnalysisTaskJetResponseV2::GetNInputTracks FAILED");
- return -1;
- }
+ // number of tracks in the event, obtained via jet finder
+
+ Int_t nInputTracks = 0;
+
+ TString jbname(fJetBranchName[1]);
+ //needs complete event, use jets without background subtraction
+ for(Int_t i=1; i<=3; ++i){
+ if(jbname.Contains(Form("B%d",i))) jbname.ReplaceAll(Form("B%d",i),"B0");
+ }
+ // use only HI event
+ if(jbname.Contains("AODextraonly")) jbname.ReplaceAll("AODextraonly","AOD");
+ else if(jbname.Contains("AODMCextraonly")) jbname.ReplaceAll("AODextraonly","AOD");
+ if(jbname.Contains("AODextra")) jbname.ReplaceAll("AODextra","AOD");
+ else if(jbname.Contains("AODMCextra")) jbname.ReplaceAll("AODextra","AOD");
+
+ if(fDebug) Printf("Multiplicity from jet branch %s", jbname.Data());
+ TClonesArray *tmpAODjets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(jbname.Data()));
+ if(!tmpAODjets && fAODOut) tmpAODjets = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(jbname.Data()));
+ if(!tmpAODjets && fAODExtension) tmpAODjets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(jbname.Data()));
+ if(!tmpAODjets){
+ Printf("Jet branch %s not found", jbname.Data());
+ Printf("AliAnalysisTaskJetResponseV2::GetNInputTracks FAILED");
+ return -1;
+ }
- for (Int_t iJet=0; iJet<tmpAODjets->GetEntriesFast(); iJet++){
- AliAODJet *jet = dynamic_cast<AliAODJet*>((*tmpAODjets)[iJet]);
- if(!jet) continue;
- TRefArray *trackList = jet->GetRefTracks();
- Int_t nTracks = trackList->GetEntriesFast();
- nInputTracks += nTracks;
- if(fDebug) Printf("#jet%d: %d tracks", iJet, nTracks);
- }
- if(fDebug) Printf("---> input tracks: %d", nInputTracks);
-
- return nInputTracks;
+ for (Int_t iJet=0; iJet<tmpAODjets->GetEntriesFast(); iJet++){
+ AliAODJet *jet = dynamic_cast<AliAODJet*>((*tmpAODjets)[iJet]);
+ if(!jet) continue;
+ TRefArray *trackList = jet->GetRefTracks();
+ Int_t nTracks = trackList->GetEntriesFast();
+ nInputTracks += nTracks;
+ if(fDebug) Printf("#jet%d: %d tracks", iJet, nTracks);
+ }
+ if(fDebug) Printf("---> input tracks: %d", nInputTracks);
+
+ return nInputTracks;
}
THnSparse* AliAnalysisTaskJetResponseV2::NewTHnSparseF(const char* name, UInt_t entries, UInt_t opt)
{
- Int_t count = 0;
- UInt_t tmp = entries;
- while(tmp!=0){
- count++;
- tmp = tmp &~ -tmp; // clear lowest bit
- }
-
- TString hnTitle(name);
- const Int_t dim = count;
- Int_t nbins[dim];
- Double_t xmin[dim];
- Double_t xmax[dim];
-
- Int_t i=0;
- Int_t c=0;
- while(c<dim && i<32){
- if(entries&(1<<i)){
- Bool_t highres = opt&(1<<i);
- TString label("");
- GetDimParams(i, highres, label, nbins[c], xmin[c], xmax[c]);
- hnTitle += Form(";%s",label.Data());
- c++;
- }
+ // generate new THnSparseF, axes are defined in GetDimParams()
+
+ Int_t count = 0;
+ UInt_t tmp = entries;
+ while(tmp!=0){
+ count++;
+ tmp = tmp &~ -tmp; // clear lowest bit
+ }
+
+ TString hnTitle(name);
+ const Int_t dim = count;
+ Int_t nbins[dim];
+ Double_t xmin[dim];
+ Double_t xmax[dim];
+
+ Int_t i=0;
+ Int_t c=0;
+ while(c<dim && i<32){
+ if(entries&(1<<i)){
+ Bool_t highres = opt&(1<<i);
+ TString label("");
+ GetDimParams(i, highres, label, nbins[c], xmin[c], xmax[c]);
+ hnTitle += Form(";%s",label.Data());
+ c++;
+ }
- i++;
- }
- hnTitle += ";";
+ i++;
+ }
+ hnTitle += ";";
- return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
+ return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
}
void AliAnalysisTaskJetResponseV2::GetDimParams(Int_t iEntry, Bool_t hr, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
{
+ // stores label and binning of axis for THnSparse
- const Double_t pi = TMath::Pi();
+ const Double_t pi = TMath::Pi();
- switch(iEntry){
+ switch(iEntry){
- case 0:
- label = "V0 centrality (%)";
- if(hr){
- nbins = 100;
- xmin = 0.;
- xmax = 100.;
- } else {
- nbins = 8;
- xmin = 0.;
- xmax = 80.;
- }
- break;
+ case 0:
+ label = "V0 centrality (%)";
+ if(hr){
+ nbins = 100;
+ xmin = 0.;
+ xmax = 100.;
+ } else {
+ nbins = 8;
+ xmin = 0.;
+ xmax = 80.;
+ }
+ break;
- case 1:
- label = "nb. of input tracks";
- if(fIsPbPb){
- if(hr){
- nbins = 400;
- xmin = 0.;
- xmax = 4000.;
- } else {
- nbins = 40;
- xmin = 0.;
- xmax = 4000.;
- }
+ case 1:
+ label = "nb. of input tracks";
+ if(fIsPbPb){
+ if(hr){
+ nbins = 400;
+ xmin = 0.;
+ xmax = 4000.;
} else {
- nbins = 40;
- xmin = 0.;
- xmax = 400.;
+ nbins = 40;
+ xmin = 0.;
+ xmax = 4000.;
}
- break;
+ } else {
+ nbins = 40;
+ xmin = 0.;
+ xmax = 400.;
+ }
+ break;
- case 2:
- label = "event plane #psi";
- if(hr){
- nbins = 30;
- xmin = 0.;
- xmax = pi;
- } else {
- nbins = 30;
- xmin = 0.;
- xmax = pi;
- }
- break;
+ case 2:
+ label = "event plane #psi";
+ if(hr){
+ nbins = 30;
+ xmin = 0.;
+ xmax = pi;
+ } else {
+ nbins = 30;
+ xmin = 0.;
+ xmax = pi;
+ }
+ break;
- case 3:
- label = "event plane bin";
- nbins = 3;
- xmin = -.5;
- xmax = 2.5;
- break;
+ case 3:
+ label = "event plane bin";
+ nbins = 3;
+ xmin = -.5;
+ xmax = 2.5;
+ break;
- case 4:
- case 5:
- if(iEntry==4)label = "#Delta#phi(RP-jet) (probe)";
- if(iEntry==5)label = "#Delta#phi(RP-jet) (rec)";
- nbins = 48;
- xmin = -pi;
- xmax = pi;
- break;
+ case 4:
+ case 5:
+ if(iEntry==4)label = "#Delta#phi(RP-jet) (probe)";
+ if(iEntry==5)label = "#Delta#phi(RP-jet) (rec)";
+ nbins = 48;
+ xmin = -pi;
+ xmax = pi;
+ break;
- case 6:
- case 7:
- if(iEntry==6)label = "probe p_{T} (GeV/c)";
- if(iEntry==7)label = "rec p_{T} (GeV/c)";
- if(hr){
- nbins = 300;
- xmin = -50.;
- xmax = 250.;
- } else {
- nbins = 50;
- xmin = 0.;
- xmax = 250.;
- }
- break;
+ case 6:
+ case 7:
+ case 27:
+ if(iEntry==6)label = "probe p_{T} (GeV/c)";
+ if(iEntry==7 || iEntry==27)label = "rec p_{T} (GeV/c)";
+ if(hr){
+ nbins = 300;
+ xmin = -50.;
+ xmax = 250.;
+ } else {
+ nbins = 50;
+ xmin = 0.;
+ xmax = 250.;
+ }
+ break;
- case 8:
- case 9:
- if(iEntry==8)label = "probe #eta";
- if(iEntry==9)label = "rec #eta";
- if(hr){
- nbins = 56;
- xmin = -.7;
- xmax = .7;
- } else {
- nbins = 28;
- xmin = -.7;
- xmax = .7;
- }
- break;
+ case 8:
+ case 9:
+ if(iEntry==8)label = "probe #eta";
+ if(iEntry==9)label = "rec #eta";
+ if(hr){
+ nbins = 56;
+ xmin = -.7;
+ xmax = .7;
+ } else {
+ nbins = 28;
+ xmin = -.7;
+ xmax = .7;
+ }
+ break;
- case 10:
- case 11:
- if(iEntry==10)label = "probe #phi";
- if(iEntry==11)label = "rec #phi";
- if(hr){
- nbins = 90; // modulo 18 (sectors)
- xmin = 0.;
- xmax = 2*pi;
- } else {
- nbins = 25;
- xmin = 0.;
- xmax = 2*pi;
- }
- break;
+ case 10:
+ case 11:
+ if(iEntry==10)label = "probe #phi";
+ if(iEntry==11)label = "rec #phi";
+ if(hr){
+ nbins = 90; // modulo 18 (sectors)
+ xmin = 0.;
+ xmax = 2*pi;
+ } else {
+ nbins = 25;
+ xmin = 0.;
+ xmax = 2*pi;
+ }
+ break;
- case 12:
- case 13:
- if(iEntry==12)label = "probe area";
- if(iEntry==13)label = "rec area";
- if(hr){
- nbins = 100;
- xmin = 0.;
- xmax = 1.;
- } else {
- nbins = 25;
- xmin = 0.;
- xmax = 1.;
- }
- break;
+ case 12:
+ case 13:
+ case 29:
+ if(iEntry==12)label = "probe area";
+ if(iEntry==13 || iEntry==29)label = "rec area";
+ if(hr){
+ nbins = 100;
+ xmin = 0.;
+ xmax = 1.;
+ } else {
+ nbins = 25;
+ xmin = 0.;
+ xmax = 1.;
+ }
+ break;
- case 14:
- label = "#Delta p_{T}";
- if(hr){
- nbins = 241;
- xmin = -120.5;
- xmax = 120.5;
- } else {
- nbins = 101;
- xmin = -101.;
- xmax = 101.;
- }
- break;
+ case 14:
+ case 28:
+ if(iEntry==14) label = "#delta p_{T}";
+ if(iEntry==28) label = "#delta";
+ if(hr){
+ nbins = 241;
+ xmin = -120.5;
+ xmax = 120.5;
+ } else {
+ nbins = 101;
+ xmin = -101.;
+ xmax = 101.;
+ }
+ break;
- case 15:
- label = "#Delta#eta";
- if(hr){
- nbins = 51;
- xmin = -1.02;
- xmax = 1.02;
- } else {
- nbins = 51;
- xmin = -1.02;
- xmax = 1.02;
- }
- break;
+ case 15:
+ label = "#Delta#eta";
+ if(hr){
+ nbins = 51;
+ xmin = -1.02;
+ xmax = 1.02;
+ } else {
+ nbins = 51;
+ xmin = -1.02;
+ xmax = 1.02;
+ }
+ break;
- case 16:
- label = "#Delta#phi";
- if(hr){
- nbins = 45;
- xmin = -pi;
- xmax = pi;
- } else {
- nbins = 45;
- xmin = -pi;
- xmax = pi;
- }
- break;
+ case 16:
+ label = "#Delta#phi";
+ if(hr){
+ nbins = 45;
+ xmin = -pi;
+ xmax = pi;
+ } else {
+ nbins = 45;
+ xmin = -pi;
+ xmax = pi;
+ }
+ break;
- case 17:
- label = "#DeltaR";
- if(hr){
- nbins = 50;
- xmin = 0.;
- xmax = 1.;
- } else {
- nbins = 50;
- xmin = 0.;
- xmax = 1.;
- }
- break;
+ case 17:
+ label = "#DeltaR";
+ if(hr){
+ nbins = 50;
+ xmin = 0.;
+ xmax = 1.;
+ } else {
+ nbins = 50;
+ xmin = 0.;
+ xmax = 1.;
+ }
+ break;
- case 18:
- label = "#Deltaarea";
- if(hr){
- nbins = 81;
- xmin = -.81;
- xmax = .81;
- } else {
- nbins = 33;
- xmin = -.825;
- xmax = .825;
- }
- break;
+ case 18:
+ label = "#Deltaarea";
+ if(hr){
+ nbins = 81;
+ xmin = -.81;
+ xmax = .81;
+ } else {
+ nbins = 33;
+ xmin = -.825;
+ xmax = .825;
+ }
+ break;
- case 19:
- label = "fraction";
- if(hr){
- nbins = 52;
- xmin = 0.;
- xmax = 1.04;
- } else {
- nbins = 52;
- xmin = 0.;
- xmax = 1.04;
- }
- break;
+ case 19:
+ case 30:
+ label = "fraction";
+ if(hr){
+ nbins = 52;
+ xmin = 0.;
+ xmax = 1.04;
+ } else {
+ nbins = 52;
+ xmin = 0.;
+ xmax = 1.04;
+ }
+ break;
- case 20:
- label = "distance to closest rec jet";
- if(hr){
- nbins = 51;
- xmin = -0.02;
- xmax = 1.;
- } else {
- nbins = 51;
- xmin = -0.02;
- xmax = 1.;
- }
- break;
+ case 20:
+ label = "distance to closest rec jet";
+ if(hr){
+ nbins = 51;
+ xmin = -0.02;
+ xmax = 1.;
+ } else {
+ nbins = 51;
+ xmin = -0.02;
+ xmax = 1.;
+ }
+ break;
- case 21:
- label = "p_{T} of closest rec jet";
- nbins = 100;
- xmin = 0.;
- xmax = 200.;
- break;
+ case 21:
+ label = "p_{T} of closest rec jet";
+ nbins = 100;
+ xmin = 0.;
+ xmax = 200.;
+ break;
- case 22:
- label = "#rho";
- nbins = 125;
- xmin = 0.;
- xmax = 250.;
- break;
+ case 22:
+ label = "#rho";
+ nbins = 125;
+ xmin = 0.;
+ xmax = 250.;
+ break;
- case 23:
- label = "abs. correction of #rho for RP";
- nbins = 51;
- xmin = -51.;
- xmax = 51.;
- break;
+ case 23:
+ label = "abs. correction of #rho for RP";
+ nbins = 51;
+ xmin = -51.;
+ xmax = 51.;
+ break;
- case 24:
- label = "local #rho";
- nbins = 125;
- xmin = 0.;
- xmax = 250.;
- break;
+ case 24:
+ label = "local #rho";
+ nbins = 125;
+ xmin = 0.;
+ xmax = 250.;
+ break;
- case 25:
- label = "local #rho / #rho";
- nbins = 500;
- xmin = 0.;
- xmax = 5.;
- break;
+ case 25:
+ label = "local #rho / #rho";
+ nbins = 500;
+ xmin = 0.;
+ xmax = 5.;
+ break;
- case 26:
- label = "p_{T,hard} bin";
- nbins = 10;
- xmin = -.5;
- xmax = 9.5;
- break;
+ case 26:
+ label = "p_{T,hard} bin";
+ nbins = 10;
+ xmin = -.5;
+ xmax = 9.5;
+ break;
- }
+ }
}
//____________________________________________________________________________
Int_t AliAnalysisTaskJetResponseV2::GetPtHardBin(Double_t ptHard){
- const Int_t nBins = 10;
- Double_t binLimits[nBins] = { 5., 11., 21., 36., 57., 84., 117., 156., 200., 249. }; // lower limits
+ // returns pt hard bin (for LHC10e14, LHC11a1x, LHC11a2x simulations)
+
+ const Int_t nBins = 10;
+ Double_t binLimits[nBins] = { 5., 11., 21., 36., 57., 84., 117., 156., 200., 249. }; // lower limits
- Int_t bin = -1;
- while(bin<nBins-1 && binLimits[bin+1]<ptHard){
- bin++;
- }
+ Int_t bin = -1;
+ while(bin<nBins-1 && binLimits[bin+1]<ptHard){
+ bin++;
+ }
- return bin;
+ return bin;
}
//____________________________________________________________________________
Double_t AliAnalysisTaskJetResponseV2::GetPt(AliAODJet *j, Int_t mode=0){
- Double_t pt = 0.;
+ // returns jet pt, also negative pt after background subtraction if available
+
+ Double_t pt = 0.;
- if(fKeepJets && mode==1){ // background subtracted pt, can be negative
- pt = j->GetPtSubtracted(0);
- }
- else{
- pt = j->Pt(); // if negative, pt=0.01
- }
+ if(fKeepJets && mode>0){ // background subtracted pt, can be negative
+ pt = j->GetPtSubtracted(0);
+ }
+ else{
+ pt = j->Pt(); // if negative, pt=0.01
+ }
- return pt;
+ return pt;
}