#include <TMath.h>
#include <TF1.h>
#include <TF2.h>
+#include <TH3.h>
#include <TH1F.h>
#include <TH2F.h>
+#include <TH3F.h>
#include <TProfile.h>
+#include <TFile.h>
+#include <TTree.h>
+#include <TKey.h>
+#include <TSystem.h>
// aliroot includes
#include <AliAnalysisTask.h>
#include <AliAnalysisManager.h>
ClassImp(AliAnalysisTaskJetMatching)
AliAnalysisTaskJetMatching::AliAnalysisTaskJetMatching() : AliAnalysisTaskEmcalJet("AliAnalysisTaskJetMatching", kTRUE),
- fDebug(0), 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), fHistNoConstSourceJet(0), fHistNoConstTargetJet(0), fHistNoConstMatchJet(0), fProfFracPtMatched(0), fProfFracPtJets(0), fProfFracNoMatched(0), fProfFracNoJets(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) {
+ 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) {
// default constructor
ClearMatchedJetsCache();
}
//_____________________________________________________________________________
AliAnalysisTaskJetMatching::AliAnalysisTaskJetMatching(const char* name) : AliAnalysisTaskEmcalJet(name, kTRUE),
- fDebug(0), 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), fHistNoConstSourceJet(0), fHistNoConstTargetJet(0), fHistNoConstMatchJet(0), fProfFracPtMatched(0), fProfFracPtJets(0), fProfFracNoMatched(0), fProfFracNoJets(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) {
+ 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) {
// constructor
ClearMatchedJetsCache();
DefineInput(0, TChain::Class());
void AliAnalysisTaskJetMatching::ExecOnce()
{
// initialize the anaysis
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
// get the stand alone jets from the input event
fSourceJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fSourceJetsName.Data()));
if(!fSourceJets) AliFatal(Form("%s: Container with name %s not found. Aborting", GetName(), fSourceJetsName.Data()));
void AliAnalysisTaskJetMatching::UserCreateOutputObjects()
{
// create output objects
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
- // add the matched jets array to the event
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
+ // add the matched jets array to the event
fMatchedJets = new TClonesArray("AliEmcalJet");
fMatchedJets->SetName(fMatchedJetsName);
fOutputList = new TList();
fHistTargetJetPt = BookTH1F("fHistTargetJetPt", "p_{t} [GeV/c]", 150, -50, 250);
fHistMatchedJetPt = (fDoDetectorResponse) ? BookTH1F("fHistDetectorLevelJet", "p_{t}^{rec}", 150, -50, 250) : BookTH1F("fHistMatchedJetPt", "p_{t} [GeV/c]", 150, -50, 250);
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);
+ if(fDoDetectorResponse) {
+ fHistDetectorResponseProb = BookTH2F("fHistDetectorResponseProb", "(p_{t}^{det} - p_{t}^{part})/p_{t}^{part}", "p_{t}^{part}", 100, -1.5, 1., 20, 0, 200);
+ fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
+ fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
+ fOutputList->Add(fh1Xsec);
+ fh1Trials = new TH1F("fh1Trials","trials root file",1,0,1);
+ fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
+ fOutputList->Add(fh1Trials);
+ fh1AvgTrials = new TH1F("fh1AvgTrials","avg trials root file",1,0,1);
+ fh1AvgTrials->GetXaxis()->SetBinLabel(1,"#sum{avg ntrials}");
+ fOutputList->Add(fh1AvgTrials);
+ }
fHistNoConstSourceJet = BookTH1F("fHistNoConstSourceJet", "number of constituents source jets", 50, 0, 50);
fHistNoConstTargetJet = BookTH1F("fHistNoConstTargetJet", "number of constituents target jets", 50, 0, 50);
fHistNoConstMatchJet = BookTH1F("fHistNoConstMatchJet", "number of constituents matched jets", 50, 0, 50);
fOutputList->Add(fProfFracPtJets);
fProfFracNoJets = new TProfile("fProfFracNoJets", "source constituents / target constituents", 15, -50, 250);
fOutputList->Add(fProfFracNoJets);
+ switch (fMatchingScheme) {
+ case kDiJet : {
+ 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);
+ 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);
+ 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);
+ 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);
+ } break;
+ default : break;
+ }
fOutputList->Sort();
PostData(1, fOutputList);
}
TH1F* AliAnalysisTaskJetMatching::BookTH1F(const char* name, const char* x, Int_t bins, Double_t min, Double_t max, Bool_t append)
{
// book a TH1F and connect it to the output container
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
- if(!fOutputList) return 0x0;
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
+ if(append && !fOutputList) return 0x0;
TString title(name);
title += Form(";%s;[counts]", x);
TH1F* histogram = new TH1F(name, title.Data(), bins, min, max);
TH2F* 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)
{
// book a TH2F and connect it to the output container
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
- if(!fOutputList) return 0x0;
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
+ if(append && !fOutputList) return 0x0;
TString title(name);
title += Form(";%s;%s", x, y);
TH2F* histogram = new TH2F(name, title.Data(), binsx, minx, maxx, binsy, miny, maxy);
return histogram;
}
//_____________________________________________________________________________
+TH3F* 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)
+{
+ // book a TH2F and connect it to the output container
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
+ if(append && !fOutputList) return 0x0;
+ TString title(name);
+ title += Form(";%s;%s;%s", x, y, z);
+ TH3F* histogram = new TH3F(name, title.Data(), binsx, minx, maxx, binsy, miny, maxy, binsz, minz, maxz);
+ histogram->Sumw2();
+ if(append) fOutputList->Add(histogram);
+ return histogram;
+}
+//_____________________________________________________________________________
+Bool_t AliAnalysisTaskJetMatching::Notify()
+{
+ // for each file get the number of trials and pythia cross section
+ // see $ALICE_ROOT/PWGJE/AliAnalysisTaskJetProperties.cxx
+ // $ALICE_ROOT/PWG/Tools/AliAnalysisHelperJetTasks.cxx
+ // this function is implenented here temporarily to avoid introducing a dependency
+ // later on this could just be a call to a static helper function
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
+ Float_t xsection(0), ftrials(1);
+ fAvgTrials = ftrials;
+ TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
+ if(tree) {
+ TFile *curfile = tree->GetCurrentFile();
+ if (!curfile) return kFALSE;
+ TString file(curfile->GetName());
+ 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),"");
+ else file.ReplaceAll(gSystem->BaseName(file.Data()),"");
+ TFile *fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec.root"));
+ if(!fxsec) {
+ fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
+ if(fxsec) {
+ TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
+ if(key) {
+ TList *list = dynamic_cast<TList*>(key->ReadObj());
+ if(list) {
+ xsection = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
+ ftrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
+ }
+ }
+ fxsec->Close();
+ }
+ } else {
+ TTree *xtree = (TTree*)fxsec->Get("Xsection");
+ if(xtree) {
+ UInt_t _ftrials = 0;
+ Double_t _xsection = 0;
+ xtree->SetBranchAddress("xsection",&_xsection);
+ xtree->SetBranchAddress("ntrials",&_ftrials);
+ xtree->GetEntry(0);
+ ftrials = _ftrials;
+ xsection = _xsection;
+ }
+ fxsec->Close();
+ }
+ fh1Xsec->Fill("<#sigma>",xsection);
+ Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
+ if(ftrials >= nEntries && nEntries > 0.) fAvgTrials = ftrials/nEntries;
+ fh1Trials->Fill("#sum{ntrials}",ftrials);
+ }
+ return kTRUE;
+}
+//_____________________________________________________________________________
Bool_t AliAnalysisTaskJetMatching::Run()
{
// execute once for each event
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
if(!(InputEvent() && fSourceJets && fTargetJets && IsEventSelected())) return kFALSE;
+ if(fh1AvgTrials) fh1AvgTrials->Fill("#sum{avg ntrials}", fAvgTrials);
// step one: do geometric matching
switch (fMatchingScheme) {
+ // FIXME having separate dedicated functions is not necessary and historical
+ // cluttered code - should be cleaned up and merged into one
case kGeoEtaPhi : {
DoGeometricMatchingEtaPhi();
+ // break if no jet was matched
+ if(!fMatchedJetContainer[1][0]) return kTRUE;
} break;
case kGeoR : {
DoGeometricMatchingR();
+ // break if no jet was matched
+ if(!fMatchedJetContainer[1][0]) return kTRUE;
+ } break;
+ case kDiJet : {
+ DoDiJetMatching();
} break;
default : break;
}
- // break if no jet was matched
- if(!fMatchedJetContainer[1][0]) return kTRUE;
// optional step two: get a bijection (avoid duplicate matches)
if(fGetBijection) GetBijection();
// optional step three: match constituents within matched jets
// stream data to output
PostMatchedJets();
FillMatchedJetHistograms();
- if(fDebug > 0) PrintInfo();
+ #ifdef DEBUGTASK
+ PrintInfo();
+ #endif
PostData(1, fOutputList);
return kTRUE;
}
void AliAnalysisTaskJetMatching::DoGeometricMatchingEtaPhi()
{
// do geometric matching based on eta phi distance between jet centers
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
fNoMatchedJets = 0; // reset the matched jet counter
Int_t iSource(fSourceJets->GetEntriesFast()), iTarget(fTargetJets->GetEntriesFast());
for(Int_t i(0); i < iSource; i++) {
if(TMath::Abs(sourcePhi-targetPhi) < fMatchPhi) { // accept the jets as matching
Bool_t isBestMatch(kTRUE);
if(fGetBijection) { // match has been found, for bijection only store it if there's no better match
- if(fDebug > 0) printf(" > Entering first bbijection test \n");
+ #ifdef DEBUGTASK
+ printf(" > Entering first bijection test \n");
+ #endif
for(Int_t k(i); k < iSource; k++) {
AliEmcalJet* candidateSourceJet(static_cast<AliEmcalJet*>(fSourceJets->At(k)));
if(PassesCuts(candidateSourceJet, 0)) continue;
- if(fDebug > 0) printf("source distance %.2f \t candidate distance %.2f \n", GetR(sourceJet, targetJet),GetR(candidateSourceJet, targetJet));
+ #ifdef DEBUGTASK
+ printf("source distance %.2f \t candidate distance %.2f \n", GetR(sourceJet, targetJet),GetR(candidateSourceJet, targetJet));
+ #endif
if(GetR(sourceJet, targetJet) > GetR(candidateSourceJet, targetJet)) {
isBestMatch = kFALSE;
break;
}
}
- if(fDebug > 0) (isBestMatch) ? printf(" kept source \n ") : printf(" we can do better (rejected source) \n");
+ #ifdef DEBUGTASK
+ (isBestMatch) ? printf(" kept source \n ") : printf(" we can do better (rejected source) \n");
+ #endif
}
if(isBestMatch) {
fMatchedJetContainer[fNoMatchedJets][0] = sourceJet;
void AliAnalysisTaskJetMatching::DoGeometricMatchingR()
{
// do geometric matching based on shortest path between jet centers
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
fNoMatchedJets = 0; // reset the matched jet counter
Int_t iSource(fSourceJets->GetEntriesFast()), iTarget(fTargetJets->GetEntriesFast());
for(Int_t i(0); i < iSource; i++) {
if(GetR(sourceJet, targetJet) <= fMatchR) {
Bool_t isBestMatch(kTRUE);
if(fGetBijection) { // match has been found, for bijection only store it if there's no better match
- if(fDebug > 0) printf(" > Entering first bijection test \n");
+ #ifdef DEBUGTASK
+ printf(" > Entering first bijection test \n");
+ #endif
for(Int_t k(i); k < iSource; k++) {
AliEmcalJet* candidateSourceJet(static_cast<AliEmcalJet*>(fSourceJets->At(k)));
if(!PassesCuts(candidateSourceJet, 0)) continue;
- if(fDebug > 0) printf("source distance %.2f \t candidate distance %.2f \n", GetR(sourceJet, targetJet),GetR(candidateSourceJet, targetJet));
+ #ifdef DEBUGTASK
+ printf("source distance %.2f \t candidate distance %.2f \n", GetR(sourceJet, targetJet),GetR(candidateSourceJet, targetJet));
+ #endif
if(GetR(sourceJet, targetJet) > GetR(candidateSourceJet, targetJet)) {
isBestMatch = kFALSE;
break;
}
}
- if(fDebug > 0) (isBestMatch) ? printf(" kept source \n ") : printf(" we can do better (rejected source) \n");
+ #ifdef DEBUGTASK
+ (isBestMatch) ? printf(" kept source \n ") : printf(" we can do better (rejected source) \n");
+ #endif
}
if(isBestMatch) {
fMatchedJetContainer[fNoMatchedJets][0] = sourceJet;
}
}
//_____________________________________________________________________________
+void AliAnalysisTaskJetMatching::DoDiJetMatching()
+{
+ // match dijets. this is in a sense a 'special' mode of the task as both jet supplied jet arrays
+ // (target and source) are the same jet collection, matching will be performed on distribution in
+ // azimuth
+ // no ouptut array is produced, dedicated dijet histo's are filled in this function instead
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
+ // get the leading jet in a given acceptance (TPC is default)
+ Int_t leadingJetIndex(-1), subLeadingJetIndex(-1);
+ // retrieve the leading jet, leadingJetIndex points to the leading jet
+ AliEmcalJet* leadingJet(GetLeadingJet(fSourceJets, leadingJetIndex));
+ if(!leadingJet) return;
+ // fill phi and eta of leading jet (should always be in selected acceptance)
+ fHistDiJetLeadingJet->Fill(leadingJet->Phi(), leadingJet->Eta(), leadingJet->Pt());
+ Double_t sourcePhi(leadingJet->Phi()), targetPhi(-1);
+ // get the sub-leading jet - faster when leading jet is also provided
+ AliEmcalJet* subLeadingJet(GetSubLeadingJet(fSourceJets, leadingJetIndex, subLeadingJetIndex));
+ if(!subLeadingJet) return;
+ else { // check if the sub-leading jet is actually the away-side jet
+ targetPhi = subLeadingJet->Phi() + TMath::Pi();
+ // rotate jets to common phase
+ if(TMath::Abs(sourcePhi) > TMath::Abs(sourcePhi+TMath::TwoPi())) sourcePhi+=TMath::TwoPi();
+ if(TMath::Abs(sourcePhi) > TMath::Abs(sourcePhi-TMath::TwoPi())) sourcePhi-=TMath::TwoPi();
+ if(TMath::Abs(targetPhi) > TMath::Abs(targetPhi+TMath::TwoPi())) targetPhi+=TMath::TwoPi();
+ if(TMath::Abs(targetPhi) > TMath::Abs(targetPhi-TMath::TwoPi())) targetPhi-=TMath::TwoPi();
+ if(TMath::Abs(sourcePhi-targetPhi) < fMatchPhi) {
+ fHistDiJet->Fill(subLeadingJet->Phi(), subLeadingJet->Eta(), leadingJet->Pt());
+ fHistDiJetDPhi->Fill(sourcePhi-targetPhi, leadingJet->Pt());
+ fHistDiJetDPt->Fill(leadingJet->Pt() - subLeadingJet->Pt(), leadingJet->Pt());
+ }
+ }
+}
+//_____________________________________________________________________________
void AliAnalysisTaskJetMatching::DoConstituentMatching()
{
// match constituents within matched jets
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
if(!fTracks) {
AliFatal(Form("%s Fatal error! To do deep matching, supply jet constituents ! \n", GetName()));
return; // coverity ...
}
}
if((float)overlap/(float)iSJ < fMinFracRecoveredConstituents || targetPt/sourceJet->Pt() < fMinFracRecoveredConstituentPt) {
- if(fDebug > 0) 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);
+ #ifdef DEBUGTASK
+ 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);
+ #endif
fMatchedJetContainer[i][0] = 0x0;
fMatchedJetContainer[i][1] = 0x0;
continue;
fProfFracNoMatched->Fill(sourceJet->Pt()-sourceRho, (double)overlap / (double)sourceJet->GetNumberOfTracks());
fProfFracNoJets->Fill(sourceJet->Pt()-sourceRho, (double)targetJet->GetNumberOfTracks() / (double)sourceJet->GetNumberOfTracks());
}
- if(fDebug > 0) {
- printf("\n > Jet A: %i const\t", iSJ);
- printf(" > Jet B %i const\t", iTJ);
- printf(" -> OVERLAP: %i tracks <- \n", overlap);
- }
+ #ifdef DEBUGTASK
+ if(fDebug > 0) {
+ printf("\n > Jet A: %i const\t", iSJ);
+ printf(" > Jet B %i const\t", iTJ);
+ printf(" -> OVERLAP: %i tracks <- \n", overlap);
}
+ #endif
+ }
}
}
//_____________________________________________________________________________
void AliAnalysisTaskJetMatching::GetBijection()
{
// bijection of source and matched jets, based on closest distance between jets
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
for(Int_t i(0); i < fNoMatchedJets; i++) {
for(Int_t j(i+1); j < fNoMatchedJets; j++) {
if(fMatchedJetContainer[i][0] == fMatchedJetContainer[j][0]) {
fMatchedJetContainer[j][0] = 0x0; // clear duplicate jet from cache
fMatchedJetContainer[j][1] = 0x0;
}
- if(fDebug > 0) printf(" found duplicate jet, chose %.2f over %.2f \n" , (rB > rA) ? rA : rB, (rB > rA) ? rB : rA);
+ #ifdef DEBUGTASK
+ printf(" found duplicate jet, chose %.2f over %.2f \n" , (rB > rA) ? rA : rB, (rB > rA) ? rB : rA);
+ #endif
}
}
}
void AliAnalysisTaskJetMatching::PostMatchedJets()
{
// post matched jets
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
fMatchedJets->Delete(); // should be a NULL operation, but added just in case
for(Int_t i(0), p(0); i < fNoMatchedJets; i++) {
if(fMatchedJetContainer[i][1]) { // duplicate jets will have NULL value here and are skipped
void AliAnalysisTaskJetMatching::FillAnalysisSummaryHistogram() const
{
// fill the analysis summary histrogram, saves all relevant analysis settigns
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
fHistAnalysisSummary->GetXaxis()->SetBinLabel(1, "fUseScaledRho");
fHistAnalysisSummary->SetBinContent(1, (int)fUseScaledRho);
fHistAnalysisSummary->GetXaxis()->SetBinLabel(2, "fMatchingScheme");
void AliAnalysisTaskJetMatching::FillMatchedJetHistograms()
{
// fill matched jet histograms
- if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #ifdef DEBUGTASK
+ printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+ #endif
for(Int_t i(0); i < fSourceJets->GetEntriesFast(); i++) {
AliEmcalJet* source = static_cast<AliEmcalJet*>(fSourceJets->At(i));
if(!source) continue;
if(fDoDetectorResponse) {
fProfFracPtJets->Fill(fMatchedJetContainer[i][0]->Pt()-sourceRho, (fMatchedJetContainer[i][1]->Pt()-targetRho) / (fMatchedJetContainer[i][0]->Pt()-sourceRho));
fProfFracNoJets->Fill(fMatchedJetContainer[i][0]->Pt()-sourceRho, (double)fMatchedJetContainer[i][1]->GetNumberOfTracks() / (double)fMatchedJetContainer[i][0]->GetNumberOfTracks());
+ fHistDetectorResponseProb->Fill((fMatchedJetContainer[i][1]->Pt()-fMatchedJetContainer[i][0]->Pt())/fMatchedJetContainer[i][0]->Pt(), fMatchedJetContainer[i][0]->Pt());
}
}
}
}
//_____________________________________________________________________________
+AliEmcalJet* AliAnalysisTaskJetMatching::GetLeadingJet(TClonesArray* source, Int_t &leadingJetIndex, Double_t etaMin, Double_t etaMax)
+{
+ // return the leading jet within giiven acceptance
+ Int_t iJets(source->GetEntriesFast());
+ Double_t pt(0);
+ AliEmcalJet* leadingJet(0x0);
+ for(Int_t i(0); i < iJets; i++) {
+ AliEmcalJet* jet = static_cast<AliEmcalJet*>(source->At(i));
+ if(jet->Eta() < etaMin || jet->Eta() > etaMax) continue;
+ if(jet->Pt() > pt) {
+ leadingJet = jet;
+ pt = leadingJet->Pt();
+ leadingJetIndex = i;
+ }
+ }
+ return leadingJet;
+}
+//_____________________________________________________________________________
+AliEmcalJet* AliAnalysisTaskJetMatching::GetSubLeadingJet(TClonesArray* source, Int_t leadingJetIndex, Int_t &subLeadingJetIndex)
+{
+ // return the sub-leading jet within given acceptance
+ // same as GetLeadingJet() but skips the leading jet (so returned jet is
+ // sub-leading by design)
+ // there is no eta requirement on the location of the sub-leading jet
+ Int_t iJets(source->GetEntriesFast());
+ // if the leading jet isn't given, retrieve it
+ if(leadingJetIndex < 0) GetLeadingJet(source, leadingJetIndex);
+ AliEmcalJet *leadingJet(0x0), *prevLeadingJet(static_cast<AliEmcalJet*>(source->At(leadingJetIndex)));
+ if(!prevLeadingJet) return 0x0;
+ Double_t pt(0), leadingPt(prevLeadingJet->Pt());
+ for(Int_t i(0); i < iJets; i++) {
+ if(i == leadingJetIndex) continue;
+ AliEmcalJet* jet = static_cast<AliEmcalJet*>(source->At(i));
+ // check if jet is actually sub-leading
+ if(jet->Pt() > pt && jet->Pt() <= leadingPt) {
+ leadingJet = jet;
+ pt = leadingJet->Pt();
+ subLeadingJetIndex = i;
+ }
+ }
+ return leadingJet;
+}
+//_____________________________________________________________________________
void AliAnalysisTaskJetMatching::PrintInfo() const
{
// print some info
printf("\n > No. of source jets from %s \n \t %i \n ", fSourceJetsName.Data(), fSourceJets->GetEntriesFast());
printf(" > No. of target jets from %s \n \t %i \n ", fTargetJetsName.Data(), fTargetJets->GetEntriesFast());
printf(" > No. of matched jets from %s \n \t %i \n ", fMatchedJetsName.Data(), fMatchedJets->GetEntriesFast());
- if(fDebug > 3) InputEvent()->GetList()->ls();
}
//_____________________________________________________________________________
void AliAnalysisTaskJetMatching::Terminate(Option_t *)
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff