changes submitted by user saiola
authormcosenti <mcosenti@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 14 Jun 2012 15:43:51 +0000 (15:43 +0000)
committermcosenti <mcosenti@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 14 Jun 2012 15:43:51 +0000 (15:43 +0000)
18 files changed:
PWGGA/EMCALJetTasks/AliAnalysisTaskEmcalJet.cxx
PWGGA/EMCALJetTasks/AliAnalysisTaskEmcalJet.h
PWGGA/EMCALJetTasks/AliAnalysisTaskRho.cxx
PWGGA/EMCALJetTasks/AliAnalysisTaskSAJF.cxx
PWGGA/EMCALJetTasks/AliAnalysisTaskSAJF.h
PWGGA/EMCALJetTasks/AliAnalysisTaskSAQA.cxx
PWGGA/EMCALJetTasks/AliAnalysisTaskScale.cxx
PWGGA/EMCALJetTasks/AliAnalysisTaskScale.h
PWGGA/EMCALJetTasks/AliEmcalJet.h
PWGGA/EMCALJetTasks/AliEmcalJetTask.cxx
PWGGA/EMCALJetTasks/AliJetEmbeddingTask.cxx
PWGGA/EMCALJetTasks/AliJetModelBaseTask.cxx
PWGGA/EMCALJetTasks/AliJetModelBaseTask.h
PWGGA/EMCALJetTasks/AliJetResponseMaker.cxx
PWGGA/EMCALJetTasks/macros/AddTaskScale.C
PWGGA/EMCALTasks/AliAnalysisTaskEmcal.cxx
PWGGA/EMCALTasks/AliAnalysisTaskEmcal.h
PWGGA/EMCALTasks/AliEmcalPicoTrackMaker.cxx

index 2a4421d..b035506 100644 (file)
@@ -102,22 +102,22 @@ void AliAnalysisTaskEmcalJet::Init()
 {
   // Init the analysis.
 
-  AliAnalysisTaskEmcal::Init();
-
   if (fAnaType == kTPC) {
     SetEtaLimits(-0.9 + fJetRadius, 0.9 - fJetRadius);
     SetPhiLimits(-10, 10);
   } else if (fAnaType == kEMCAL) {
     AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance();
-    if (!geom) {
-      AliFatal(Form("%s: Can not create geometry", GetName()));
-      return;
+    if (geom) {
+      SetEtaLimits(geom->GetArm1EtaMin() + fJetRadius, geom->GetArm1EtaMax() - fJetRadius);
+      SetPhiLimits(geom->GetArm1PhiMin() * TMath::DegToRad() + fJetRadius, geom->GetArm1PhiMax() * TMath::DegToRad() - fJetRadius);
+    }
+    else {
+      AliWarning(Form("%s: Can not create geometry", GetName()));
     }
-    SetEtaLimits(geom->GetArm1EtaMin() + fJetRadius, geom->GetArm1EtaMax() - fJetRadius);
-    SetPhiLimits(geom->GetArm1PhiMin() * TMath::DegToRad() + fJetRadius, geom->GetArm1PhiMax() * TMath::DegToRad() - fJetRadius);
   } else {
     AliWarning(Form("%s: Analysis type not recognized! Assuming kTPC!", GetName()));
     SetAnaType(kTPC);
+    Init();
   }
 
   SetInitialized();
@@ -164,7 +164,13 @@ Bool_t AliAnalysisTaskEmcalJet::RetrieveEventObjects()
   if (!fJetsName.IsNull() && !fJets) {
     fJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fJetsName));
     if (!fJets) {
-      AliWarning(Form("%s: Could not retrieve jets %s!", GetName(), fJetsName.Data())); 
+      AliError(Form("%s: Could not retrieve jets %s!", GetName(), fJetsName.Data()));
+      return kFALSE;
+    }
+    else if (!fJets->GetClass()->GetBaseClass("AliEmcalJet")) {
+      AliError(Form("%s: Collection %s does not contain AliEmcalJet objects!", GetName(), fJetsName.Data())); 
+      fJets = 0;
+      return kFALSE;
     }
   }
 
index ff4ce0e..4f59cfa 100644 (file)
@@ -29,14 +29,14 @@ class AliAnalysisTaskEmcalJet : public AliAnalysisTaskEmcal {
   void                        SetPhiLimits(Float_t min, Float_t max)               { fMinPhi = min, fMaxPhi = max ; }
   void                        SetPtBiasJetClus(Float_t b)                          { fPtBiasJetClus  = b          ; }
   void                        SetPtBiasJetTrack(Float_t b)                         { fPtBiasJetTrack = b          ; }
+  void                        Init();
  
  protected:
+  Bool_t                      RetrieveEventObjects();
   Bool_t                      AcceptJet(AliEmcalJet* jet, Bool_t bias = kTRUE, Bool_t upCut = kTRUE)   const;
   Bool_t                      AcceptBiasJet(AliEmcalJet* jet)                                          const;
-  void                        Init();
   Bool_t                      IsJetTrack(AliEmcalJet* jet, Int_t itrack, Bool_t sorted = kTRUE)        const;
   Bool_t                      IsJetCluster(AliEmcalJet* jet, Int_t iclus, Bool_t sorted = kTRUE)       const;
-  Bool_t                      RetrieveEventObjects();
 
   Float_t                     fJetRadius;                  // jet radius
   TString                     fJetsName;                   // name of jet collection
index c268f80..74fd124 100644 (file)
@@ -164,12 +164,15 @@ void AliAnalysisTaskRho::UserCreateOutputObjects()
 void AliAnalysisTaskRho::UserExec(Option_t *) 
 {
   // Main loop, called for each event.
-  
+
   if (!fIsInit) {
     ExecOnce();
     fIsInit = 1;
   }
 
+  if (!fJets)
+    return;
+
   fRho->SetVal(-1);
   if (fRhoScaled)
     fRhoScaled->SetVal(-1);
index 983f043..0a31937 100644 (file)
@@ -21,6 +21,7 @@
 #include "AliVTrack.h"
 #include "AliEmcalJet.h"
 #include "AliVEventHandler.h"
+#include "AliRhoParameter.h"
 #include "AliLog.h"
 
 #include "AliAnalysisTaskSAJF.h"
@@ -32,16 +33,18 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() :
   AliAnalysisTaskEmcalJet("AliAnalysisTaskSAJF", kTRUE),
   fMinRC2LJ(1.0),
   fEmbJetsName("EmbJets"),
+  fEmbTracksName(""),
+  fEmbCaloName(""),
   fRandTracksName("TracksRandomized"),
   fRandCaloName("CaloClustersRandomized"),
   fRhoName("Rho"),
-  fSkipEventsWithNoJets(kTRUE),
   fEmbJets(0),
+  fEmbTracks(0),
+  fEmbCaloClusters(0),
   fRandTracks(0),
   fRandCaloClusters(0),
   fRho(0),
   fHistCentrality(0),
-  fHistRejectedEvents(0),
   fHistRhoVSleadJetPt(0),
   fHistRCPhiEta(0),
   fHistRCPtExLJVSDPhiLJ(0),
@@ -54,20 +57,25 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() :
   // Default constructor.
 
   for (Int_t i = 0; i < 4; i++) {
+    fHistEvents[i] = 0;
+    fHistTracksPt[i] = 0;
+    fHistClustersPt[i] = 0;
     fHistJetPhiEta[i] = 0;
     fHistJetsPt[i] = 0;
     fHistJetsPtArea[i] = 0;
-    fHistJetsPtNonBias[i] = 0;
-    fHistJetsPtAreaNonBias[i] = 0;
     fHistLeadingJetPt[i] = 0;
     fHist2LeadingJetPt[i] = 0;
     fHistJetsNEFvsPt[i] = 0;
     fHistJetsZvsPt[i] = 0;
+    fHistMaxTrackPtvsJetPt[i] = 0;
+    fHistMaxClusPtvsJetPt[i] = 0;
+    fHistMaxPartPtvsJetPt[i] = 0;
+    fHistMaxTrackPtvsJetCorrPt[i] = 0;
+    fHistMaxClusPtvsJetCorrPt[i] = 0;
+    fHistMaxPartPtvsJetCorrPt[i] = 0;
     fHistRho[i] = 0;
     fHistCorrJetsPt[i] = 0;
     fHistCorrJetsPtArea[i] = 0;
-    fHistCorrJetsPtNonBias[i] = 0;
-    fHistCorrJetsPtAreaNonBias[i] = 0;
     fHistCorrLeadingJetPt[i] = 0;
     fHistRCPt[i] = 0;
     fHistRCPtExLJ[i] = 0;
@@ -75,7 +83,8 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() :
     fHistDeltaPtRC[i] = 0;
     fHistDeltaPtRCExLJ[i] = 0;
     fHistDeltaPtRCRand[i] = 0;
-    fHistEmbJets[i] = 0;
+    fHistEmbJetsPt[i] = 0;
+    fHistEmbJetsCorrPt[i] = 0;
     fHistEmbPart[i] = 0;
     fHistDeltaPtEmb[i] = 0;
   }
@@ -86,16 +95,18 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) :
   AliAnalysisTaskEmcalJet(name, kTRUE),
   fMinRC2LJ(1.0),
   fEmbJetsName("EmbJets"),
+  fEmbTracksName(""),
+  fEmbCaloName(""),
   fRandTracksName("TracksRandomized"),
   fRandCaloName("CaloClustersRandomized"),
   fRhoName("Rho"),
-  fSkipEventsWithNoJets(kTRUE),
   fEmbJets(0),
+  fEmbTracks(0),
+  fEmbCaloClusters(0),
   fRandTracks(0),
   fRandCaloClusters(0),
   fRho(0),
   fHistCentrality(0),
-  fHistRejectedEvents(0),
   fHistRhoVSleadJetPt(0),
   fHistRCPhiEta(0),
   fHistRCPtExLJVSDPhiLJ(0),
@@ -107,20 +118,25 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) :
   // Standard constructor.
 
   for (Int_t i = 0; i < 4; i++) {
+    fHistEvents[i] = 0;
+    fHistTracksPt[i] = 0;
+    fHistClustersPt[i] = 0;
     fHistJetPhiEta[i] = 0;
     fHistJetsPt[i] = 0;
     fHistJetsPtArea[i] = 0;
-    fHistJetsPtNonBias[i] = 0;
-    fHistJetsPtAreaNonBias[i] = 0;
     fHistLeadingJetPt[i] = 0;
     fHist2LeadingJetPt[i] = 0;
     fHistJetsNEFvsPt[i] = 0;
     fHistJetsZvsPt[i] = 0;
+    fHistMaxTrackPtvsJetPt[i] = 0;
+    fHistMaxClusPtvsJetPt[i] = 0;
+    fHistMaxPartPtvsJetPt[i] = 0;
+    fHistMaxTrackPtvsJetCorrPt[i] = 0;
+    fHistMaxClusPtvsJetCorrPt[i] = 0;
+    fHistMaxPartPtvsJetCorrPt[i] = 0;
     fHistRho[i] = 0;
     fHistCorrJetsPt[i] = 0;
     fHistCorrJetsPtArea[i] = 0;
-    fHistCorrJetsPtNonBias[i] = 0;
-    fHistCorrJetsPtAreaNonBias[i] = 0;
     fHistCorrLeadingJetPt[i] = 0;
     fHistRCPt[i] = 0;
     fHistRCPtExLJ[i] = 0;
@@ -128,7 +144,8 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) :
     fHistDeltaPtRC[i] = 0;
     fHistDeltaPtRCExLJ[i] = 0;
     fHistDeltaPtRCRand[i] = 0;
-    fHistEmbJets[i] = 0;
+    fHistEmbJetsPt[i] = 0;
+    fHistEmbJetsCorrPt[i] = 0;
     fHistEmbPart[i] = 0;
     fHistDeltaPtEmb[i] = 0;
   }
@@ -158,15 +175,6 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
   fHistCentrality->GetYaxis()->SetTitle("counts");
   fOutput->Add(fHistCentrality);
 
-  fHistRejectedEvents = new TH1F("fHistRejectedEvents","fHistRejectedEvents", 6, 0, 6);
-  fHistRejectedEvents->GetXaxis()->SetTitle("Reason");
-  fHistRejectedEvents->GetYaxis()->SetTitle("counts");
-  fHistRejectedEvents->GetXaxis()->SetBinLabel(1, "Rho <= 0");
-  fHistRejectedEvents->GetXaxis()->SetBinLabel(2, "Max Jet <= 0");
-  fHistRejectedEvents->GetXaxis()->SetBinLabel(3, "Max Jet not found");
-  fHistRejectedEvents->GetXaxis()->SetBinLabel(4, "No jets");
-  fOutput->Add(fHistRejectedEvents);
-
   fHistRhoVSleadJetPt = new TH2F("fHistRhoVSleadJetPt","fHistRhoVSleadJetPt", fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
   fHistRhoVSleadJetPt->GetXaxis()->SetTitle("#rho * area [GeV/c]");
   fHistRhoVSleadJetPt->GetYaxis()->SetTitle("Leading jet p_{T} [GeV/c]");
@@ -187,24 +195,54 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
   fHistRhoVSRCPt->GetYaxis()->SetTitle("rigid cone p_{T} [GeV/c]");
   fOutput->Add(fHistRhoVSRCPt);
 
-  fHistEmbJetPhiEta = new TH2F("fHistEmbJetPhiEta","Phi-Eta distribution of embedded jets", 40, -2, 2, 64, 0, 6.4);
-  fHistEmbJetPhiEta->GetXaxis()->SetTitle("#eta");
-  fHistEmbJetPhiEta->GetYaxis()->SetTitle("#phi");
-  fOutput->Add(fHistEmbJetPhiEta);
-
-  fHistEmbPartPhiEta = new TH2F("fHistEmbPartPhiEta","Phi-Eta distribution of embedded particles", 40, -2, 2, 64, 0, 6.4);
-  fHistEmbPartPhiEta->GetXaxis()->SetTitle("#eta");
-  fHistEmbPartPhiEta->GetYaxis()->SetTitle("#phi");
-  fOutput->Add(fHistEmbPartPhiEta);
-
-  fHistRhoVSEmbBkg = new TH2F("fHistRhoVSEmbBkg","fHistRhoVSEmbBkg", fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
-  fHistRhoVSEmbBkg->GetXaxis()->SetTitle("rho * area [GeV/c]");
-  fHistRhoVSEmbBkg->GetYaxis()->SetTitle("background of embedded track [GeV/c]");
-  fOutput->Add(fHistRhoVSEmbBkg);
+  if (!fEmbJetsName.IsNull()) {
+    fHistEmbJetPhiEta = new TH2F("fHistEmbJetPhiEta","Phi-Eta distribution of embedded jets", 40, -2, 2, 64, 0, 6.4);
+    fHistEmbJetPhiEta->GetXaxis()->SetTitle("#eta");
+    fHistEmbJetPhiEta->GetYaxis()->SetTitle("#phi");
+    fOutput->Add(fHistEmbJetPhiEta);
+    
+    fHistEmbPartPhiEta = new TH2F("fHistEmbPartPhiEta","Phi-Eta distribution of embedded particles", 40, -2, 2, 64, 0, 6.4);
+    fHistEmbPartPhiEta->GetXaxis()->SetTitle("#eta");
+    fHistEmbPartPhiEta->GetYaxis()->SetTitle("#phi");
+    fOutput->Add(fHistEmbPartPhiEta);
+    
+    fHistRhoVSEmbBkg = new TH2F("fHistRhoVSEmbBkg","fHistRhoVSEmbBkg", fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
+    fHistRhoVSEmbBkg->GetXaxis()->SetTitle("rho * area [GeV/c]");
+    fHistRhoVSEmbBkg->GetYaxis()->SetTitle("background of embedded track [GeV/c]");
+    fOutput->Add(fHistRhoVSEmbBkg);
+  }
 
   TString histname;
 
   for (Int_t i = 0; i < 4; i++) {
+    histname = "fHistEvents_";
+    histname += i;
+    fHistEvents[i] = new TH1F(histname,histname, 6, 0, 6);
+    fHistEvents[i]->GetXaxis()->SetTitle("Event state");
+    fHistEvents[i]->GetYaxis()->SetTitle("counts");
+    fHistEvents[i]->GetXaxis()->SetBinLabel(1, "Rho <= 0");
+    fHistEvents[i]->GetXaxis()->SetBinLabel(2, "No jets");
+    fHistEvents[i]->GetXaxis()->SetBinLabel(3, "Max Jet not found");
+    fHistEvents[i]->GetXaxis()->SetBinLabel(4, "Max Jet <= 0");
+    fHistEvents[i]->GetXaxis()->SetBinLabel(5, "Accepted");
+    fOutput->Add(fHistEvents[i]);
+
+    histname = "fHistTracksPt_";
+    histname += i;
+    fHistTracksPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins / 2.5, fMinBinPt, fMaxBinPt / 2.5);
+    fHistTracksPt[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+    fHistTracksPt[i]->GetYaxis()->SetTitle("counts");
+    fOutput->Add(fHistTracksPt[i]);
+
+    if (fAnaType == kEMCAL) {    
+      histname = "fHistClustersPt_";
+      histname += i;
+      fHistClustersPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins / 2.5, fMinBinPt, fMaxBinPt / 2.5);
+      fHistClustersPt[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+      fHistClustersPt[i]->GetYaxis()->SetTitle("counts");
+      fOutput->Add(fHistClustersPt[i]);
+    }
+
     histname = "fHistJetPhiEta_";
     histname += i;
     fHistJetPhiEta[i] = new TH2F(histname.Data(), histname.Data(), 40, -2, 2, 64, 0, 6.4);
@@ -215,42 +253,28 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
     histname = "fHistJetsPt_";
     histname += i;
     fHistJetsPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
-    fHistJetsPt[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+    fHistJetsPt[i]->GetXaxis()->SetTitle("p_{T}^{raw} [GeV/c]");
     fHistJetsPt[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistJetsPt[i]);
 
     histname = "fHistJetsPtArea_";
     histname += i;
     fHistJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
-    fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+    fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T}^{raw} [GeV/c]");
     fHistJetsPtArea[i]->GetYaxis()->SetTitle("area");
     fOutput->Add(fHistJetsPtArea[i]);
 
-    histname = "fHistJetsPtNonBias_";
-    histname += i;
-    fHistJetsPtNonBias[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
-    fHistJetsPtNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
-    fHistJetsPtNonBias[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistJetsPtNonBias[i]);
-
-    histname = "fHistJetsPtAreaNonBias_";
-    histname += i;
-    fHistJetsPtAreaNonBias[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
-    fHistJetsPtAreaNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
-    fHistJetsPtAreaNonBias[i]->GetYaxis()->SetTitle("area");
-    fOutput->Add(fHistJetsPtAreaNonBias[i]);
-
     histname = "fHistLeadingJetPt_";
     histname += i;
     fHistLeadingJetPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
-    fHistLeadingJetPt[i]->GetXaxis()->SetTitle("p_{T} [GeV]");
+    fHistLeadingJetPt[i]->GetXaxis()->SetTitle("p_{T}^{raw} [GeV]");
     fHistLeadingJetPt[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistLeadingJetPt[i]);
 
     histname = "fHist2LeadingJetPt_";
     histname += i;
     fHist2LeadingJetPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
-    fHist2LeadingJetPt[i]->GetXaxis()->SetTitle("p_{T} [GeV]");
+    fHist2LeadingJetPt[i]->GetXaxis()->SetTitle("p_{T}^{raw} [GeV]");
     fHist2LeadingJetPt[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHist2LeadingJetPt[i]);
 
@@ -258,7 +282,7 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
     histname += i;
     fHistJetsZvsPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
     fHistJetsZvsPt[i]->GetXaxis()->SetTitle("Z");
-    fHistJetsZvsPt[i]->GetYaxis()->SetTitle("p_{T} [GeV/c]");
+    fHistJetsZvsPt[i]->GetYaxis()->SetTitle("p_{T}^{raw} [GeV/c]");
     fOutput->Add(fHistJetsZvsPt[i]);
 
     if (fAnaType == kEMCAL) {
@@ -266,10 +290,56 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
       histname += i;
       fHistJetsNEFvsPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
       fHistJetsNEFvsPt[i]->GetXaxis()->SetTitle("NEF");
-      fHistJetsNEFvsPt[i]->GetYaxis()->SetTitle("p_{T} [GeV/c]");
+      fHistJetsNEFvsPt[i]->GetYaxis()->SetTitle("p_{T}^{raw} [GeV/c]");
       fOutput->Add(fHistJetsNEFvsPt[i]);
     }
 
+    histname = "fHistMaxTrackPtvsJetPt_";
+    histname += i;
+    fHistMaxTrackPtvsJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins / 2.5, fMinBinPt, fMaxBinPt / 2.5);
+    fHistMaxTrackPtvsJetPt[i]->GetXaxis()->SetTitle("p_{T}^{jet} [GeV/c]");
+    fHistMaxTrackPtvsJetPt[i]->GetYaxis()->SetTitle("p_{T}^{track} [GeV/c]");
+    fOutput->Add(fHistMaxTrackPtvsJetPt[i]);
+
+    if (fAnaType == kEMCAL) {
+      histname = "fHistMaxClusPtvsJetPt_";
+      histname += i;
+      fHistMaxClusPtvsJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins / 2.5, fMinBinPt, fMaxBinPt / 2.5);
+      fHistMaxClusPtvsJetPt[i]->GetXaxis()->SetTitle("p_{T}^{jet} [GeV/c]");
+      fHistMaxClusPtvsJetPt[i]->GetYaxis()->SetTitle("p_{T}^{clus} [GeV/c]");
+      fOutput->Add(fHistMaxClusPtvsJetPt[i]);
+    }
+
+    histname = "fHistMaxPartPtvsJetPt_";
+    histname += i;
+    fHistMaxPartPtvsJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins / 2.5, fMinBinPt, fMaxBinPt / 2.5);
+    fHistMaxPartPtvsJetPt[i]->GetXaxis()->SetTitle("p_{T}^{jet} [GeV/c]");
+    fHistMaxPartPtvsJetPt[i]->GetYaxis()->SetTitle("p_{T}^{part} [GeV/c]");
+    fOutput->Add(fHistMaxPartPtvsJetPt[i]);
+
+    histname = "fHistMaxTrackPtvsJetCorrPt_";
+    histname += i;
+    fHistMaxTrackPtvsJetCorrPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt, fNbins / 2.5, fMinBinPt, fMaxBinPt / 2.5);
+    fHistMaxTrackPtvsJetCorrPt[i]->GetXaxis()->SetTitle("p_{T}^{jet} [GeV/c]");
+    fHistMaxTrackPtvsJetCorrPt[i]->GetYaxis()->SetTitle("p_{T}^{track} [GeV/c]");
+    fOutput->Add(fHistMaxTrackPtvsJetCorrPt[i]);
+
+    if (fAnaType == kEMCAL) {
+      histname = "fHistMaxClusPtvsJetCorrPt_";
+      histname += i;
+      fHistMaxClusPtvsJetCorrPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt, fNbins / 2.5, fMinBinPt, fMaxBinPt / 2.5);
+      fHistMaxClusPtvsJetCorrPt[i]->GetXaxis()->SetTitle("p_{T}^{jet} [GeV/c]");
+      fHistMaxClusPtvsJetCorrPt[i]->GetYaxis()->SetTitle("p_{T}^{clus} [GeV/c]");
+      fOutput->Add(fHistMaxClusPtvsJetCorrPt[i]);
+    }
+
+    histname = "fHistMaxPartPtvsJetCorrPt_";
+    histname += i;
+    fHistMaxPartPtvsJetCorrPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt, fNbins / 2.5, fMinBinPt, fMaxBinPt / 2.5);
+    fHistMaxPartPtvsJetCorrPt[i]->GetXaxis()->SetTitle("p_{T}^{jet} [GeV/c]");
+    fHistMaxPartPtvsJetCorrPt[i]->GetYaxis()->SetTitle("p_{T}^{part} [GeV/c]");
+    fOutput->Add(fHistMaxPartPtvsJetCorrPt[i]);
+
     histname = "fHistRho_";
     histname += i;
     fHistRho[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
@@ -280,35 +350,21 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
     histname = "fHistCorrJetsPt_";
     histname += i;
     fHistCorrJetsPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt);
-    fHistCorrJetsPt[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+    fHistCorrJetsPt[i]->GetXaxis()->SetTitle("p_{T}^{corr} [GeV/c]");
     fHistCorrJetsPt[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistCorrJetsPt[i]);
 
     histname = "fHistCorrJetsPtArea_";
     histname += i;
     fHistCorrJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
-    fHistCorrJetsPtArea[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+    fHistCorrJetsPtArea[i]->GetXaxis()->SetTitle("p_{T}^{corr} [GeV/c]");
     fHistCorrJetsPtArea[i]->GetYaxis()->SetTitle("area");
     fOutput->Add(fHistCorrJetsPtArea[i]);
 
-    histname = "fHistCorrJetsPtNonBias_";
-    histname += i;
-    fHistCorrJetsPtNonBias[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt);
-    fHistCorrJetsPtNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
-    fHistCorrJetsPtNonBias[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistCorrJetsPtNonBias[i]);
-
-    histname = "fHistCorrJetsPtAreaNonBias_";
-    histname += i;
-    fHistCorrJetsPtAreaNonBias[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
-    fHistCorrJetsPtAreaNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
-    fHistCorrJetsPtAreaNonBias[i]->GetYaxis()->SetTitle("area");
-    fOutput->Add(fHistCorrJetsPtAreaNonBias[i]);
-
     histname = "fHistCorrLeadingJetPt_";
     histname += i;
     fHistCorrLeadingJetPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt);
-    fHistCorrLeadingJetPt[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+    fHistCorrLeadingJetPt[i]->GetXaxis()->SetTitle("p_{T}^{RC} [GeV/c]");
     fHistCorrLeadingJetPt[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistCorrLeadingJetPt[i]);
     
@@ -322,58 +378,67 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
     histname = "fHistRCPtExLJ_";
     histname += i;
     fHistRCPtExLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
-    fHistRCPtExLJ[i]->GetXaxis()->SetTitle("rigid cone p_{T} [GeV/c]");
+    fHistRCPtExLJ[i]->GetXaxis()->SetTitle("rigid cone p_{T}^{RC} [GeV/c]");
     fHistRCPtExLJ[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistRCPtExLJ[i]);
 
     histname = "fHistRCPtRand_";
     histname += i;
     fHistRCPtRand[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
-    fHistRCPtRand[i]->GetXaxis()->SetTitle("rigid cone p_{T} [GeV/c]");
+    fHistRCPtRand[i]->GetXaxis()->SetTitle("rigid cone p_{T}^{RC} [GeV/c]");
     fHistRCPtRand[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistRCPtRand[i]);
 
     histname = "fHistDeltaPtRC_";
     histname += i;
     fHistDeltaPtRC[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt - fMaxBinPt / 2 + binWidth / 2, fMinBinPt + fMaxBinPt / 2 + binWidth / 2);
-    fHistDeltaPtRC[i]->GetXaxis()->SetTitle("#deltap_{T} [GeV/c]");
+    fHistDeltaPtRC[i]->GetXaxis()->SetTitle("#deltap_{T}^{RC} [GeV/c]");
     fHistDeltaPtRC[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistDeltaPtRC[i]);
 
     histname = "fHistDeltaPtRCExLJ_";
     histname += i;
     fHistDeltaPtRCExLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt - fMaxBinPt / 2 + binWidth / 2, fMinBinPt + fMaxBinPt / 2 + binWidth / 2);
-    fHistDeltaPtRCExLJ[i]->GetXaxis()->SetTitle("#deltap_{T} [GeV/c]");
+    fHistDeltaPtRCExLJ[i]->GetXaxis()->SetTitle("#deltap_{T}^{RC} [GeV/c]");
     fHistDeltaPtRCExLJ[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistDeltaPtRCExLJ[i]);
 
     histname = "fHistDeltaPtRCRand_";
     histname += i;
     fHistDeltaPtRCRand[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt - fMaxBinPt / 2 + binWidth / 2, fMinBinPt + fMaxBinPt / 2 + binWidth / 2);
-    fHistDeltaPtRCRand[i]->GetXaxis()->SetTitle("#deltap_{T} [GeV/c]");
+    fHistDeltaPtRCRand[i]->GetXaxis()->SetTitle("#deltap_{T}^{RC} [GeV/c]");
     fHistDeltaPtRCRand[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistDeltaPtRCRand[i]);
 
-    histname = "fHistEmbJets_";
-    histname += i;
-    fHistEmbJets[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
-    fHistEmbJets[i]->GetXaxis()->SetTitle("embedded jet p_{T} [GeV/c]");
-    fHistEmbJets[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistEmbJets[i]);
-
-    histname = "fHistEmbPart_";
-    histname += i;
-    fHistEmbPart[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
-    fHistEmbPart[i]->GetXaxis()->SetTitle("embedded particle p_{T} [GeV/c]");
-    fHistEmbPart[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistEmbPart[i]);
+    if (!fEmbJetsName.IsNull()) {
+      histname = "fHistEmbJetsPt_";
+      histname += i;
+      fHistEmbJetsPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
+      fHistEmbJetsPt[i]->GetXaxis()->SetTitle("embedded jet p_{T}^{raw} [GeV/c]");
+      fHistEmbJetsPt[i]->GetYaxis()->SetTitle("counts");
+      fOutput->Add(fHistEmbJetsPt[i]);
 
-    histname = "fHistDeltaPtEmb_";
-    histname += i;
-    fHistDeltaPtEmb[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt - fMaxBinPt / 2 + binWidth / 2, fMinBinPt + fMaxBinPt / 2 + binWidth / 2);
-    fHistDeltaPtEmb[i]->GetXaxis()->SetTitle("#deltap_{T} [GeV/c]");
-    fHistDeltaPtEmb[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistDeltaPtEmb[i]);
+      histname = "fHistEmbJetsCorrPt_";
+      histname += i;
+      fHistEmbJetsCorrPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
+      fHistEmbJetsCorrPt[i]->GetXaxis()->SetTitle("embedded jet p_{T}^{corr} [GeV/c]");
+      fHistEmbJetsCorrPt[i]->GetYaxis()->SetTitle("counts");
+      fOutput->Add(fHistEmbJetsCorrPt[i]);
+      
+      histname = "fHistEmbPart_";
+      histname += i;
+      fHistEmbPart[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
+      fHistEmbPart[i]->GetXaxis()->SetTitle("embedded particle p_{T}^{emb} [GeV/c]");
+      fHistEmbPart[i]->GetYaxis()->SetTitle("counts");
+      fOutput->Add(fHistEmbPart[i]);
+      
+      histname = "fHistDeltaPtEmb_";
+      histname += i;
+      fHistDeltaPtEmb[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt - fMaxBinPt / 2 + binWidth / 2, fMinBinPt + fMaxBinPt / 2 + binWidth / 2);
+      fHistDeltaPtEmb[i]->GetXaxis()->SetTitle("#deltap_{T}^{emb} [GeV/c]");
+      fHistDeltaPtEmb[i]->GetYaxis()->SetTitle("counts");
+      fOutput->Add(fHistDeltaPtEmb[i]);
+    }
   }
 
   PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
@@ -386,38 +451,77 @@ Bool_t AliAnalysisTaskSAJF::RetrieveEventObjects()
 
   if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
     return kFALSE;
-
-  fRho = -1;
   
-  if (strcmp(fRhoName,"")) {
-    TParameter<Double_t> *rho = dynamic_cast<TParameter<Double_t>*>(InputEvent()->FindListObject(fRhoName));
-  
-    if (rho) {
-      fRho = rho->GetVal();
-    }
-    else {
-      AliWarning(Form("Could not retrieve rho %s!", fRhoName.Data()));
+  if (!fRhoName.IsNull() && !fRho) {
+    fRho = dynamic_cast<AliRhoParameter*>(InputEvent()->FindListObject(fRhoName));
+    if (!fRho) {
+      AliError(Form("%s: Could not retrieve rho %s!", GetName(), fRhoName.Data()));
+      return kFALSE;
     }
   }
   
-  if (strcmp(fEmbJetsName,"")) {
+  if (!fEmbJetsName.IsNull() && !fEmbJets) {
     fEmbJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fEmbJetsName));
     if (!fEmbJets) {
-      AliWarning(Form("Could not retrieve emb jets %s!", fEmbJetsName.Data())); 
+      AliError(Form("%s: Could not retrieve emb jets %s!", GetName(), fEmbJetsName.Data()));
+      return kFALSE;
+    }
+    else if (!fEmbJets->GetClass()->GetBaseClass("AliEmcalJet")) {
+      AliError(Form("%s: Collection %s does not contain AliEmcalJet objects!", GetName(), fEmbJetsName.Data())); 
+      fEmbJets = 0;
+      return kFALSE;
     }
   }
 
-  if (strcmp(fRandCaloName,"") && fAnaType == kEMCAL) {
+  if (!fEmbCaloName.IsNull() && fAnaType == kEMCAL && !fEmbCaloClusters) {
+    fEmbCaloClusters =  dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fEmbCaloName));
+    if (!fEmbCaloClusters) {
+      AliError(Form("%s: Could not retrieve embedded clusters %s!", GetName(), fEmbCaloName.Data()));
+      return kFALSE;
+    }
+    else if (!fEmbCaloClusters->GetClass()->GetBaseClass("AliVCluster") && !fEmbCaloClusters->GetClass()->GetBaseClass("AliEmcalParticle")) {
+      AliError(Form("%s: Collection %s does not contain AliVCluster nor AliEmcalParticle objects!", GetName(), fEmbCaloName.Data())); 
+      fEmbCaloClusters = 0;
+      return kFALSE;
+    }
+  }
+
+  if (!fEmbTracksName.IsNull() && !fEmbTracks) {
+    fEmbTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fEmbTracksName));
+    if (!fEmbTracks) {
+      AliError(Form("%s: Could not retrieve embedded tracks %s!", GetName(), fEmbTracksName.Data()));
+      return kFALSE;
+    }
+    else if (!fEmbTracks->GetClass()->GetBaseClass("AliVParticle") && !fEmbTracks->GetClass()->GetBaseClass("AliEmcalParticle")) {
+      AliError(Form("%s: Collection %s does not contain AliVParticle nor AliEmcalParticle objects!", GetName(), fEmbTracksName.Data())); 
+      fEmbTracks = 0;
+      return kFALSE;
+    }
+  }
+
+  if (!fRandCaloName.IsNull() && fAnaType == kEMCAL && !fRandCaloClusters) {
     fRandCaloClusters =  dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fRandCaloName));
     if (!fRandCaloClusters) {
-      AliWarning(Form("Could not retrieve randomized clusters %s!", fRandCaloName.Data())); 
+      AliError(Form("%s: Could not retrieve randomized clusters %s!", GetName(), fRandCaloName.Data()));
+      return kFALSE;
+    }
+    else if (!fRandCaloClusters->GetClass()->GetBaseClass("AliVCluster") && !fRandCaloClusters->GetClass()->GetBaseClass("AliEmcalParticle")) {
+      AliError(Form("%s: Collection %s does not contain AliVCluster nor AliEmcalParticle objects!", GetName(), fRandCaloName.Data())); 
+      fRandCaloClusters = 0;
+      return kFALSE;
     }
   }
 
-  if (strcmp(fRandTracksName,"")) {
+  if (!fRandTracksName.IsNull() && !fRandTracks) {
     fRandTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fRandTracksName));
     if (!fRandTracks) {
-      AliWarning(Form("Could not retrieve randomized tracks %s!", fRandTracksName.Data())); 
+      AliError(Form("%s: Could not retrieve randomized tracks %s!", GetName(), fRandTracksName.Data()));
+      return kFALSE;
+    }
+    else if (!fRandTracks->GetClass()->GetBaseClass("AliVParticle") && !fRandTracks->GetClass()->GetBaseClass("AliEmcalParticle")) {
+      AliError(Form("%s: Collection %s does not contain AliVParticle nor AliEmcalParticle objects!", GetName(), fRandTracksName.Data())); 
+      fRandTracks = 0;
+      return kFALSE;
     }
   }
 
@@ -429,52 +533,59 @@ Bool_t AliAnalysisTaskSAJF::FillHistograms()
 {
   // Fill histograms.
 
-  if (fRho < 0) {
-    fHistRejectedEvents->Fill("Rho <= 0", 1);
-    return kFALSE;
+  // Before filling any histogram, check if the event is interesting
+  Double_t rho = fRho->GetVal();
+
+  if (rho < 0) {  // rho < 0, probably a diffractive event, skipping
+    fHistEvents[fCentBin]->Fill("Rho <= 0", 1);
+    return kTRUE;
   }
 
   Int_t maxJetIndex  = -1;
   Int_t max2JetIndex = -1;
 
-  // ************
-  // General histograms
-  // _________________________________
-
-  DoJetLoop();
-
   GetLeadingJets(maxJetIndex, max2JetIndex);
   
-  if (maxJetIndex < 0) {
-    fHistRejectedEvents->Fill("No jets", 1);
-    return kFALSE;
+  if (maxJetIndex < 0) { // no accept jet, skipping
+    fHistEvents[fCentBin]->Fill("No jets", 1);
+    return kTRUE;
   }
 
-  AliEmcalJet* jet = dynamic_cast<AliEmcalJet*>(fJets->At(maxJetIndex));
-  if (!jet) {
-    fHistRejectedEvents->Fill("Max Jet not found", 1);
-    return kFALSE;
+  AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(maxJetIndex));
+  if (!jet) {  // error, I cannot get the lead jet from collection (should never happen), skipping
+    fHistEvents[fCentBin]->Fill("Max Jet not found", 1);
+    return kTRUE;
   }
 
-  Float_t maxJetCorrPt = jet->Pt() - fRho * jet->Area();
+  // OK, event accepted
 
-  if (fSkipEventsWithNoJets && maxJetCorrPt <= 0) {
-    fHistRejectedEvents->Fill("Max Jet <= 0", 1);
-    return kFALSE;
-  }
+  Float_t maxJetCorrPt = jet->Pt() - rho * jet->Area();
+
+  if (maxJetCorrPt <= 0)
+    fHistEvents[fCentBin]->Fill("Max Jet <= 0", 1);
+
+  fHistEvents[fCentBin]->Fill("Accepted", 1);
+
+  // ************
+  // General histograms
+  // _________________________________
+
+  DoJetLoop();
+  DoTrackLoop();
+  DoClusterLoop();
 
   fHistCentrality->Fill(fCent);
-  fHistRho[fCentBin]->Fill(fRho);
+  fHistRho[fCentBin]->Fill(rho);
 
   if (jet) {
     fHistLeadingJetPt[fCentBin]->Fill(jet->Pt());
-    fHistRhoVSleadJetPt->Fill(fRho * jet->Area(), jet->Pt());
+    fHistRhoVSleadJetPt->Fill(rho * jet->Area(), jet->Pt());
     fHistCorrLeadingJetPt[fCentBin]->Fill(maxJetCorrPt);
   }
   
   AliEmcalJet* jet2 = 0;
   if (max2JetIndex >= 0)
-    jet2 = dynamic_cast<AliEmcalJet*>(fJets->At(max2JetIndex));
+    jet2 = static_cast<AliEmcalJet*>(fJets->At(max2JetIndex));
 
   if (jet2)
     fHist2LeadingJetPt[fCentBin]->Fill(jet2->Pt());
@@ -484,25 +595,26 @@ Bool_t AliAnalysisTaskSAJF::FillHistograms()
   // _________________________________
   
   const Float_t rcArea = fJetRadius * fJetRadius * TMath::Pi();
-
+  const Bool_t IsMCEvent = (Bool_t)(fTracksName.Contains("Randomized") || fTracksName.Contains("Embedded"));
+  
   // Simple random cones
   Float_t RCpt = 0;
   Float_t RCeta = 0;
   Float_t RCphi = 0;
-  GetRigidCone(RCpt, RCeta, RCphi, kFALSE, 0);
+  GetRigidCone(RCpt, RCeta, RCphi, IsMCEvent, 0);
   if (RCpt > 0) {
     fHistRCPt[fCentBin]->Fill(RCpt / rcArea);
-    fHistDeltaPtRC[fCentBin]->Fill(RCpt - rcArea * fRho);
+    fHistDeltaPtRC[fCentBin]->Fill(RCpt - rcArea * rho);
   }
   
   // Random cones far from leading jet
   Float_t RCptExLJ = 0;
   Float_t RCetaExLJ = 0;
   Float_t RCphiExLJ = 0;
-  GetRigidCone(RCptExLJ, RCetaExLJ, RCphiExLJ, kFALSE, jet);
+  GetRigidCone(RCptExLJ, RCetaExLJ, RCphiExLJ, IsMCEvent, jet);
   if (RCptExLJ > 0) {
     fHistRCPhiEta->Fill(RCetaExLJ, RCphiExLJ);
-    fHistRhoVSRCPt->Fill(fRho, RCptExLJ / rcArea);
+    fHistRhoVSRCPt->Fill(rho, RCptExLJ / rcArea);
 
     Float_t dphi = RCphiExLJ - jet->Phi();
     if (dphi > 4.8) dphi -= TMath::Pi() * 2;
@@ -510,7 +622,7 @@ Bool_t AliAnalysisTaskSAJF::FillHistograms()
     fHistRCPtExLJVSDPhiLJ->Fill(RCptExLJ, dphi);
     
     fHistRCPtExLJ[fCentBin]->Fill(RCptExLJ / rcArea);
-    fHistDeltaPtRCExLJ[fCentBin]->Fill(RCptExLJ - rcArea * fRho);
+    fHistDeltaPtRCExLJ[fCentBin]->Fill(RCptExLJ - rcArea * rho);
   }
 
   // Random cones with randomized particles
@@ -520,7 +632,7 @@ Bool_t AliAnalysisTaskSAJF::FillHistograms()
   GetRigidCone(RCptRand, RCetaRand, RCphiRand, kTRUE, 0, fRandTracks, fRandCaloClusters);
   if (RCptRand > 0) {
     fHistRCPtRand[fCentBin]->Fill(RCptRand / rcArea);
-    fHistDeltaPtRCRand[fCentBin]->Fill(RCptRand - rcArea * fRho);
+    fHistDeltaPtRCRand[fCentBin]->Fill(RCptRand - rcArea * rho);
   }
 
   // ************
@@ -531,19 +643,20 @@ Bool_t AliAnalysisTaskSAJF::FillHistograms()
     return kTRUE;
 
   AliEmcalJet *embJet  = 0;
-  TObject     *maxPart = 0;
+  TObject     *embPart = 0;
 
-  DoEmbJetLoop(embJet, maxPart);
+  DoEmbJetLoop(embJet, embPart);
 
   if (embJet) {
 
-    fHistEmbJets[fCentBin]->Fill(embJet->Pt());
+    fHistEmbJetsPt[fCentBin]->Fill(embJet->Pt());
+    fHistEmbJetsCorrPt[fCentBin]->Fill(embJet->Pt() - rho * embJet->Area());
     fHistEmbPart[fCentBin]->Fill(embJet->MCPt());
     fHistEmbJetPhiEta->Fill(embJet->Eta(), embJet->Phi());
-    fHistDeltaPtEmb[fCentBin]->Fill(embJet->Pt() - embJet->Area() * fRho - embJet->MCPt());
-    fHistRhoVSEmbBkg->Fill(embJet->Area() * fRho, embJet->Pt() - embJet->MCPt());
+    fHistDeltaPtEmb[fCentBin]->Fill(embJet->Pt() - embJet->Area() * rho - embJet->MCPt());
+    fHistRhoVSEmbBkg->Fill(embJet->Area() * rho, embJet->Pt() - embJet->MCPt());
 
-    AliVCluster *cluster = dynamic_cast<AliVCluster*>(maxPart);
+    AliVCluster *cluster = dynamic_cast<AliVCluster*>(embPart);
     if (cluster) {
       Float_t pos[3];
       cluster->GetPosition(pos);
@@ -552,7 +665,7 @@ Bool_t AliAnalysisTaskSAJF::FillHistograms()
       fHistEmbPartPhiEta->Fill(clusVec.Eta(), clusVec.Phi());
     }
     else {
-      AliVTrack *track = dynamic_cast<AliVTrack*>(maxPart);
+      AliVTrack *track = dynamic_cast<AliVTrack*>(embPart);
       if (track) {
        fHistEmbPartPhiEta->Fill(track->Eta(), track->Phi());
       }
@@ -578,23 +691,25 @@ void AliAnalysisTaskSAJF::GetLeadingJets(Int_t &maxJetIndex, Int_t &max2JetIndex
   if (!fJets)
     return;
 
+  const Double_t rho = fRho->GetVal();
+
   const Int_t njets = fJets->GetEntriesFast();
 
   Float_t maxJetPt = -999;
   Float_t max2JetPt = -999;
   for (Int_t ij = 0; ij < njets; ij++) {
 
-    AliEmcalJet* jet = dynamic_cast<AliEmcalJet*>(fJets->At(ij));
+    AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(ij));
 
     if (!jet) {
       AliError(Form("Could not receive jet %d", ij));
       continue;
     }  
 
-    if (!AcceptJet(jet, fSkipEventsWithNoJets))
+    if (!AcceptJet(jet))
       continue;
 
-    Float_t corrPt = jet->Pt() - fRho * jet->Area();
+    Float_t corrPt = jet->Pt() - rho * jet->Area();
 
     if (maxJetIndex == -1 || maxJetPt < corrPt) {
       max2JetPt = maxJetPt;
@@ -610,6 +725,59 @@ void AliAnalysisTaskSAJF::GetLeadingJets(Int_t &maxJetIndex, Int_t &max2JetIndex
 }
 
 //________________________________________________________________________
+void AliAnalysisTaskSAJF::DoClusterLoop()
+{
+  // Do cluster loop.
+
+  if (!fCaloClusters)
+    return;
+
+  Int_t nclusters =  fCaloClusters->GetEntriesFast();
+
+  for (Int_t iClusters = 0; iClusters < nclusters; iClusters++) {
+    AliVCluster* cluster = static_cast<AliVCluster*>(fCaloClusters->At(iClusters));
+    if (!cluster) {
+      AliError(Form("Could not receive cluster %d", iClusters));
+      continue;
+    }  
+
+    if (!AcceptCluster(cluster, kTRUE)) 
+      continue;
+
+    fHistClustersPt[fCentBin]->Fill(cluster->E());
+  }
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskSAJF::DoTrackLoop()
+{
+  // Do track loop.
+
+  if (!fTracks)
+    return;
+
+  Int_t ntracks = fTracks->GetEntriesFast();
+
+  for (Int_t i = 0; i < ntracks; i++) {
+
+    AliVParticle* track = static_cast<AliVParticle*>(fTracks->At(i)); // pointer to reconstructed to track  
+
+    if (!track) {
+      AliError(Form("Could not retrieve track %d",i)); 
+      continue; 
+    }
+
+    AliVTrack* vtrack = dynamic_cast<AliVTrack*>(track); 
+    
+    if (vtrack && !AcceptTrack(vtrack, kTRUE)) 
+      continue;
+    
+    fHistTracksPt[fCentBin]->Fill(track->Pt());
+  }
+}
+
+
+//________________________________________________________________________
 void AliAnalysisTaskSAJF::DoJetLoop()
 {
   // Do the jet loop.
@@ -617,39 +785,42 @@ void AliAnalysisTaskSAJF::DoJetLoop()
   if (!fJets)
     return;
 
+  const Double_t rho = fRho->GetVal();
+
   const Int_t njets = fJets->GetEntriesFast();
 
   for (Int_t ij = 0; ij < njets; ij++) {
 
-    AliEmcalJet* jet = dynamic_cast<AliEmcalJet*>(fJets->At(ij));
+    AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(ij));
 
     if (!jet) {
       AliError(Form("Could not receive jet %d", ij));
       continue;
     }  
 
-    if (!AcceptJet(jet, kFALSE))
+    if (!AcceptJet(jet))
       continue;
 
-    Float_t corrPt = jet->Pt() - fRho * jet->Area();
-    
-    fHistJetsPtNonBias[fCentBin]->Fill(jet->Pt());
-    fHistJetsPtAreaNonBias[fCentBin]->Fill(jet->Pt(), jet->Area());
-    fHistCorrJetsPtNonBias[fCentBin]->Fill(corrPt);
-    fHistCorrJetsPtAreaNonBias[fCentBin]->Fill(corrPt, jet->Area());
-
-    if (!AcceptBiasJet(jet))
-      continue;
+    Float_t corrPt = jet->Pt() - rho * jet->Area();
 
     fHistJetsPt[fCentBin]->Fill(jet->Pt());
     fHistJetsPtArea[fCentBin]->Fill(jet->Pt(), jet->Area());
     fHistCorrJetsPt[fCentBin]->Fill(corrPt);
     fHistCorrJetsPtArea[fCentBin]->Fill(corrPt, jet->Area());
 
+    fHistMaxTrackPtvsJetPt[fCentBin]->Fill(jet->Pt(), jet->MaxTrackPt());
+    fHistMaxPartPtvsJetPt[fCentBin]->Fill(jet->Pt(), jet->MaxPartPt());
+
+    fHistMaxTrackPtvsJetCorrPt[fCentBin]->Fill(corrPt, jet->MaxTrackPt());
+    fHistMaxPartPtvsJetCorrPt[fCentBin]->Fill(corrPt, jet->MaxPartPt());
+
     fHistJetPhiEta[fCentBin]->Fill(jet->Eta(), jet->Phi());
 
-    if (fAnaType == kEMCAL)
+    if (fAnaType == kEMCAL) {
+      fHistMaxClusPtvsJetPt[fCentBin]->Fill(jet->Pt(), jet->MaxClusterPt());
+      fHistMaxClusPtvsJetCorrPt[fCentBin]->Fill(corrPt, jet->MaxClusterPt());
       fHistJetsNEFvsPt[fCentBin]->Fill(jet->NEF(), jet->Pt());
+    }
 
     if (fTracks) {
       for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
@@ -674,7 +845,7 @@ void AliAnalysisTaskSAJF::DoJetLoop()
 }
 
 //________________________________________________________________________
-void AliAnalysisTaskSAJF::DoEmbJetLoop(AliEmcalJet* &embJet, TObject* &maxPart)
+void AliAnalysisTaskSAJF::DoEmbJetLoop(AliEmcalJet* &embJet, TObject* &embPart)
 {
   // Do the embedded jet loop.
 
@@ -687,7 +858,7 @@ void AliAnalysisTaskSAJF::DoEmbJetLoop(AliEmcalJet* &embJet, TObject* &maxPart)
 
   for (Int_t ij = 0; ij < nembjets; ij++) {
       
-    AliEmcalJet* jet = dynamic_cast<AliEmcalJet*>(fEmbJets->At(ij));
+    AliEmcalJet* jet = static_cast<AliEmcalJet*>(fEmbJets->At(ij));
       
     if (!jet) {
       AliError(Form("Could not receive jet %d", ij));
@@ -703,7 +874,7 @@ void AliAnalysisTaskSAJF::DoEmbJetLoop(AliEmcalJet* &embJet, TObject* &maxPart)
     AliVParticle *maxTrack = 0;
 
     for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
-      AliVParticle *track = jet->TrackAt(it, fTracks);
+      AliVParticle *track = jet->TrackAt(it, fEmbTracks);
       
       if (!track) 
        continue;
@@ -718,7 +889,7 @@ void AliAnalysisTaskSAJF::DoEmbJetLoop(AliEmcalJet* &embJet, TObject* &maxPart)
     AliVCluster *maxClus = 0;
 
     for (Int_t ic = 0; ic < jet->GetNumberOfClusters(); ic++) {
-      AliVCluster *cluster = jet->ClusterAt(ic, fCaloClusters);
+      AliVCluster *cluster = jet->ClusterAt(ic, fEmbCaloClusters);
       
       if (!cluster) 
        continue;
@@ -736,11 +907,11 @@ void AliAnalysisTaskSAJF::DoEmbJetLoop(AliEmcalJet* &embJet, TObject* &maxPart)
     }
 
     if ((maxClus && maxTrack && maxClusVect.Et() > maxTrack->Pt()) || (maxClus && !maxTrack)) {
-      maxPart = maxClus;
+      embPart = maxClus;
       embJet = jet;
     }
     else if (maxTrack) {
-      maxPart = maxTrack;
+      embPart = maxTrack;
       embJet = jet;
     }
 
@@ -792,21 +963,22 @@ void AliAnalysisTaskSAJF::GetRigidCone(Float_t &pt, Float_t &eta, Float_t &phi,
     repeats++;
   } while (dLJ < fMinRC2LJ && repeats < 999);
 
-  if (repeats == 999)
+  if (repeats == 999) {
+    AliWarning("Could not get random cone!");
     return;
+  }
 
   if (fAnaType == kEMCAL && clusters) {
     Int_t nclusters =  clusters->GetEntriesFast();
     for (Int_t iClusters = 0; iClusters < nclusters; iClusters++) {
-      AliVCluster* cluster = dynamic_cast<AliVCluster*>(clusters->At(iClusters));
+      AliVCluster* cluster = static_cast<AliVCluster*>(clusters->At(iClusters));
       if (!cluster) {
        AliError(Form("Could not receive cluster %d", iClusters));
        continue;
       }  
       
-      if (!(cluster->IsEMCAL())) continue;
-      
-      if (!AcceptCluster(cluster, acceptMC)) continue;
+      if (!AcceptCluster(cluster, acceptMC))
+       continue;
       
       TLorentzVector nPart;
       cluster->GetMomentum(nPart, const_cast<Double_t*>(fVertex));
@@ -825,13 +997,14 @@ void AliAnalysisTaskSAJF::GetRigidCone(Float_t &pt, Float_t &eta, Float_t &phi,
   if (tracks) {
     Int_t ntracks = tracks->GetEntriesFast();
     for(Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
-      AliVTrack* track = dynamic_cast<AliVTrack*>(tracks->At(iTracks));         
+      AliVTrack* track = static_cast<AliVTrack*>(tracks->At(iTracks));         
       if(!track) {
        AliError(Form("Could not retrieve track %d",iTracks)); 
        continue; 
       }
       
-      if (!AcceptTrack(track, acceptMC)) continue;
+      if (!AcceptTrack(track, acceptMC)) 
+       continue;
       
       Float_t tracketa = track->Eta();
       Float_t trackphi = track->Phi();
@@ -842,7 +1015,6 @@ void AliAnalysisTaskSAJF::GetRigidCone(Float_t &pt, Float_t &eta, Float_t &phi,
        trackphi -= 2 * TMath::Pi();
       
       Float_t d = TMath::Sqrt((tracketa - eta) * (tracketa - eta) + (trackphi - phi) * (trackphi - phi));
-    
       if (d <= fJetRadius)
        pt += track->Pt();
     }
@@ -855,6 +1027,17 @@ void AliAnalysisTaskSAJF::Init()
 
   AliAnalysisTaskEmcalJet::Init();
 
+  if (!fEmbJetsName.IsNull()) {
+    if (fEmbTracksName.IsNull()) {
+      fEmbTracksName = fTracksName;
+      fEmbTracksName += "Embedded";
+    }
+    if (fEmbCaloName.IsNull() && fAnaType == kEMCAL) {
+      fEmbCaloName = fCaloName;
+      fEmbCaloName += "Embedded";
+    }
+  }
+
   const Float_t semiDiag = TMath::Sqrt((fMaxPhi - fMinPhi) * (fMaxPhi - fMinPhi) + 
                                        (fMaxEta - fMinEta) * (fMaxEta - fMinEta)) / 2;
   if (fMinRC2LJ > semiDiag * 0.5) {
index f3c4b9b..5722322 100644 (file)
@@ -7,6 +7,7 @@ class TClonesArray;
 class TString;
 class TH1F;
 class TH2F;
+class AliRhoParameter;
 
 #include "AliAnalysisTaskEmcalJet.h"
 
@@ -25,8 +26,9 @@ class AliAnalysisTaskSAJF : public AliAnalysisTaskEmcalJet {
   void                        SetEmbJetsName(const char *n)                        { fEmbJetsName             = n          ; }
   void                        SetRandTracksName(const char *n)                     { fRandTracksName          = n          ; }
   void                        SetRandClusName(const char *n)                       { fRandCaloName            = n          ; }
-  void                        SetRhoName(const char *n)                            { fRhoName                 = n          ; }
-  void                        SetSkipEventsWithNoJets(Bool_t s)                    { fSkipEventsWithNoJets    = s          ; } 
+  void                        SetEmbTracksName(const char *n)                      { fEmbTracksName          = n          ; }
+  void                        SetEmbClusName(const char *n)                        { fEmbCaloName            = n          ; }
+  void                        SetRhoName(const char *n)                            { fRhoName                 = n          ; } 
 
  protected:
 
@@ -34,44 +36,51 @@ class AliAnalysisTaskSAJF : public AliAnalysisTaskEmcalJet {
   Bool_t                      FillHistograms()                                                                              ;
   void                        GetLeadingJets(Int_t &maxJetIndex, Int_t &max2JetIndex)                                       ;
   void                        DoJetLoop()                                                                                   ;
-  void                        DoEmbJetLoop(AliEmcalJet* &embJet, TObject* &maxPart)                                         ;
-  void                        DoTrackLoop(Int_t maxJetIndex)                                                                ;
-  void                        DoClusterLoop(Int_t maxJetIndex)                                                              ;
+  void                        DoEmbJetLoop(AliEmcalJet* &embJet, TObject* &embPart)                                         ;
+  void                        DoTrackLoop()                                                                                 ;
+  void                        DoClusterLoop()                                                                               ;
   void                        GetRigidCone(Float_t &pt, Float_t &eta, Float_t &phi, Bool_t acceptMC = kFALSE, 
                                           AliEmcalJet *jet = 0, TClonesArray* tracks = 0, TClonesArray* clusters = 0) const;
 
   Float_t                     fMinRC2LJ;                   // Minimum distance random cone to leading jet
-  TString                     fEmbJetsName;                // Name of embedded jets collection
+  TString                     fEmbJetsName;                // Name of embedded jet collection
+  TString                     fEmbTracksName;              // Name of embedded track collection
+  TString                     fEmbCaloName;                // Name of embedded calo cluster collection
   TString                     fRandTracksName;             // Name of randomized track collection
   TString                     fRandCaloName;               // Name of randomized calo cluster collection
   TString                     fRhoName;                    // Name of rho object
-  Bool_t                      fSkipEventsWithNoJets;       // Whether or not skip events with no jets
 
   TClonesArray               *fEmbJets;                    //!Embedded Jets
+  TClonesArray               *fEmbTracks;                  //!Embedded tracks
+  TClonesArray               *fEmbCaloClusters;            //!Embedded clusters  
   TClonesArray               *fRandTracks;                 //!Randomized tracks
   TClonesArray               *fRandCaloClusters;           //!Randomized clusters
-  Float_t                     fRho;                        //!Event rho
+  AliRhoParameter            *fRho;                        //!Event rho
 
   // General histograms
   TH1F                       *fHistCentrality;             //!Event centrality distribution
-  TH1F                       *fHistRejectedEvents;         //!Rejected events
+  TH1F                       *fHistEvents[4];              //!Events accepted/rejected
+  TH1F                       *fHistTracksPt[4];            //!Inclusive track pt spectrum
+  TH1F                       *fHistClustersPt[4];          //!Inclusive clusters pt spectrum
   TH2F                       *fHistJetPhiEta[4];           //!Phi-Eta distribution of jets
   TH1F                       *fHistJetsPt[4];              //!Inclusive jet pt spectrum
   TH2F                       *fHistJetsPtArea[4];          //!Jet pt vs. area
-  TH1F                       *fHistJetsPtNonBias[4];       //!Non-biased inclusive jet pt spectrum
-  TH2F                       *fHistJetsPtAreaNonBias[4];   //!Non-biased jet pt vs. area
   TH1F                       *fHistLeadingJetPt[4];        //!Leading jet pt spectrum
   TH1F                       *fHist2LeadingJetPt[4];       //!Second leading jet pt spectrum
   TH2F                       *fHistJetsNEFvsPt[4];         //!Jet neutral energy fraction vs. jet pt
   TH2F                       *fHistJetsZvsPt[4];           //!Constituent Pt over Jet Pt ratio vs. jet pt
+  TH2F                       *fHistMaxTrackPtvsJetPt[4];   //!Max constituent track pt vs. jet pt
+  TH2F                       *fHistMaxClusPtvsJetPt[4];    //!Max constituent cluster pt vs. jet pt
+  TH2F                       *fHistMaxPartPtvsJetPt[4];    //!Max constituent particle (track or cluster) pt vs. jet pt
+  TH2F                       *fHistMaxTrackPtvsJetCorrPt[4];   //!Max constituent track pt vs. jet pt
+  TH2F                       *fHistMaxClusPtvsJetCorrPt[4];    //!Max constituent cluster pt vs. jet pt
+  TH2F                       *fHistMaxPartPtvsJetCorrPt[4];    //!Max constituent particle (track or cluster) pt vs. jet pt
 
   // Rho
   TH1F                       *fHistRho[4];                    //!Rho distribution
   TH2F                       *fHistRhoVSleadJetPt;            //!Area(leadjetarea) * rho vs. leading jet pt
   TH1F                       *fHistCorrJetsPt[4];             //!Corrected inclusive jet pt spectrum
   TH2F                       *fHistCorrJetsPtArea[4];         //!Jet pt vs. area
-  TH1F                       *fHistCorrJetsPtNonBias[4];      //!Corrected non-biased inclusive jet pt spectrum
-  TH2F                       *fHistCorrJetsPtAreaNonBias[4];  //!Non-biased jet pt vs. area
   TH1F                       *fHistCorrLeadingJetPt[4];       //!Corrected leading jet pt spectrum
 
   // Random cones
@@ -86,7 +95,8 @@ class AliAnalysisTaskSAJF : public AliAnalysisTaskEmcalJet {
   TH1F                       *fHistDeltaPtRCRand[4];       //!deltaPt = Pt(RC) - A * rho, randomzied particles
 
   // Jet embedding
-  TH1F                       *fHistEmbJets[4];             //!Pt distribution of embedded jets
+  TH1F                       *fHistEmbJetsPt[4];           //!Pt distribution of embedded jets
+  TH1F                       *fHistEmbJetsCorrPt[4];       //!Pt distribution of embedded jets
   TH1F                       *fHistEmbPart[4];             //!Pt distribution of embedded particle
   TH2F                       *fHistEmbJetPhiEta;           //!Phi-Eta distribution of embedded jets
   TH2F                       *fHistEmbPartPhiEta;          //!Phi-Eta distribution of embedded particles
index 7a2bb35..87ea537 100644 (file)
@@ -181,11 +181,6 @@ void AliAnalysisTaskSAQA::UserCreateOutputObjects()
     fHistClusPhiEta->GetYaxis()->SetTitle("#phi");
     fOutput->Add(fHistClusPhiEta);
   }
-       
-  fHistJetsPhiEta = new TH2F("fHistJetsPhiEta","Phi-Eta distribution of jets", 80, -2, 2, 128, 0, 6.4);
-  fHistJetsPhiEta->GetXaxis()->SetTitle("#eta");
-  fHistJetsPhiEta->GetYaxis()->SetTitle("#phi");
-  fOutput->Add(fHistJetsPhiEta);
 
   if (fAnaType == kEMCAL) {
    
@@ -238,52 +233,60 @@ void AliAnalysisTaskSAQA::UserCreateOutputObjects()
   fHistTrackPhi[4]->SetLineColor(kBlack);
   fHistTrackEta[4]->SetLineColor(kBlack);
 
-  TString histname;
+  if (!fJetsName.IsNull()) {
 
-  for (Int_t i = 0; i < 4; i++) {
-    histname = "fHistJetsPtNonBias_";
-    histname += i;
-    fHistJetsPtNonBias[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5);
-    fHistJetsPtNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
-    fHistJetsPtNonBias[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistJetsPtNonBias[i]);
-
-    histname = "fHistJetsPtTrack_";
-    histname += i;
-    fHistJetsPtTrack[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5);
-    fHistJetsPtTrack[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
-    fHistJetsPtTrack[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistJetsPtTrack[i]);
-
-    if (fAnaType == kEMCAL) {
-      histname = "fHistJetsPtClus_";
+    fHistJetsPhiEta = new TH2F("fHistJetsPhiEta","Phi-Eta distribution of jets", 80, -2, 2, 128, 0, 6.4);
+    fHistJetsPhiEta->GetXaxis()->SetTitle("#eta");
+    fHistJetsPhiEta->GetYaxis()->SetTitle("#phi");
+    fOutput->Add(fHistJetsPhiEta);
+
+    TString histname;
+
+    for (Int_t i = 0; i < 4; i++) {
+      histname = "fHistJetsPtNonBias_";
       histname += i;
-      fHistJetsPtClus[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5);
-      fHistJetsPtClus[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
-      fHistJetsPtClus[i]->GetYaxis()->SetTitle("counts");
-      fOutput->Add(fHistJetsPtClus[i]);
-    }
+      fHistJetsPtNonBias[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5);
+      fHistJetsPtNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+      fHistJetsPtNonBias[i]->GetYaxis()->SetTitle("counts");
+      fOutput->Add(fHistJetsPtNonBias[i]);
+
+      histname = "fHistJetsPtTrack_";
+      histname += i;
+      fHistJetsPtTrack[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5);
+      fHistJetsPtTrack[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+      fHistJetsPtTrack[i]->GetYaxis()->SetTitle("counts");
+      fOutput->Add(fHistJetsPtTrack[i]);
+
+      if (fAnaType == kEMCAL) {
+       histname = "fHistJetsPtClus_";
+       histname += i;
+       fHistJetsPtClus[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5);
+       fHistJetsPtClus[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+       fHistJetsPtClus[i]->GetYaxis()->SetTitle("counts");
+       fOutput->Add(fHistJetsPtClus[i]);
+      }
+
+      histname = "fHistJetsPt_";
+      histname += i;
+      fHistJetsPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5);
+      fHistJetsPt[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+      fHistJetsPt[i]->GetYaxis()->SetTitle("counts");
+      fOutput->Add(fHistJetsPt[i]);
 
-    histname = "fHistJetsPt_";
-    histname += i;
-    fHistJetsPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5);
-    fHistJetsPt[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
-    fHistJetsPt[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistJetsPt[i]);
-
-    histname = "fHistJetsPtAreaNonBias_";
-    histname += i;
-    fHistJetsPtAreaNonBias[i] = new TH2F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
-    fHistJetsPtAreaNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
-    fHistJetsPtAreaNonBias[i]->GetYaxis()->SetTitle("area");
-    fOutput->Add(fHistJetsPtAreaNonBias[i]);
-
-    histname = "fHistJetsPtArea_";
-    histname += i;
-    fHistJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
-    fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
-    fHistJetsPtArea[i]->GetYaxis()->SetTitle("area");
-    fOutput->Add(fHistJetsPtArea[i]);
+      histname = "fHistJetsPtAreaNonBias_";
+      histname += i;
+      fHistJetsPtAreaNonBias[i] = new TH2F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
+      fHistJetsPtAreaNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+      fHistJetsPtAreaNonBias[i]->GetYaxis()->SetTitle("area");
+      fOutput->Add(fHistJetsPtAreaNonBias[i]);
+
+      histname = "fHistJetsPtArea_";
+      histname += i;
+      fHistJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins * 2.5, fMinBinPt, fMaxBinPt * 2.5, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
+      fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
+      fHistJetsPtArea[i]->GetYaxis()->SetTitle("area");
+      fOutput->Add(fHistJetsPtArea[i]);
+    }
   }
 
   PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
@@ -297,10 +300,19 @@ Bool_t AliAnalysisTaskSAQA::RetrieveEventObjects()
   if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
     return kFALSE;
 
-  if (strcmp(fTrgClusName,"") && fDoTrigger) {
+  if (!fTrgClusName.IsNull() && fDoTrigger && !fTrgClusters) {
     fTrgClusters =  dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTrgClusName));
     if (!fTrgClusters) {
-      AliWarning(Form("Could not retrieve trigger clusters!")); 
+      AliError(Form("%s: Could not retrieve trigger clusters %s!", GetName(), fTrgClusName.Data())); 
+      return kFALSE;
+    }
+    else {
+      TClass *cl = fTrgClusters->GetClass();
+      if (!cl->GetBaseClass("AliVCluster") && !cl->GetBaseClass("AliEmcalParticle")) {
+       AliError(Form("%s: Collection %s does not contain AliVCluster nor AliEmcalParticle objects!", GetName(), fTrgClusName.Data())); 
+       fTrgClusters = 0;
+       return kFALSE;
+      }
     }
   }
 
@@ -389,15 +401,14 @@ Float_t AliAnalysisTaskSAQA::DoClusterLoop()
   Int_t nclusters =  fCaloClusters->GetEntriesFast();
 
   for (Int_t iClusters = 0; iClusters < nclusters; iClusters++) {
-    AliVCluster* cluster = dynamic_cast<AliVCluster*>(fCaloClusters->At(iClusters));
+    AliVCluster* cluster = static_cast<AliVCluster*>(fCaloClusters->At(iClusters));
     if (!cluster) {
       AliError(Form("Could not receive cluster %d", iClusters));
       continue;
     }  
-    
-    if (!(cluster->IsEMCAL())) continue;
 
-    if (!AcceptCluster(cluster, kTRUE)) continue;
+    if (!AcceptCluster(cluster, kTRUE))
+      continue;
 
     fHistClustersEnergy->Fill(cluster->E());
 
@@ -423,7 +434,6 @@ Float_t AliAnalysisTaskSAQA::DoTrackLoop()
 
   Float_t sum = 0;
 
-  // Track loop 
   Int_t ntracks = fTracks->GetEntriesFast();
   Int_t nclusters = 0;
   if (fCaloClusters)
@@ -431,7 +441,7 @@ Float_t AliAnalysisTaskSAQA::DoTrackLoop()
 
   for (Int_t i = 0; i < ntracks; i++) {
 
-    AliVParticle* track = dynamic_cast<AliVParticle*>(fTracks->At(i)); // pointer to reconstructed to track  
+    AliVParticle* track = static_cast<AliVParticle*>(fTracks->At(i)); // pointer to reconstructed to track  
 
     if (!track) {
       AliError(Form("Could not retrieve track %d",i)); 
@@ -480,7 +490,7 @@ void AliAnalysisTaskSAQA::DoJetLoop()
 
   for (Int_t ij = 0; ij < njets; ij++) {
 
-    AliEmcalJet* jet = dynamic_cast<AliEmcalJet*>(fJets->At(ij));
+    AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(ij));
 
     if (!jet) {
       AliError(Form("Could not receive jet %d", ij));
@@ -521,7 +531,7 @@ Float_t AliAnalysisTaskSAQA::DoTriggerClusLoop()
   Float_t maxTrgClus = 0;
 
   for (Int_t iClusters = 0; iClusters < ntrgclusters; iClusters++) {
-    AliVCluster* cluster = dynamic_cast<AliVCluster*>(fTrgClusters->At(iClusters));
+    AliVCluster* cluster = static_cast<AliVCluster*>(fTrgClusters->At(iClusters));
     if (!cluster) {
       AliError(Form("Could not receive cluster %d", iClusters));
       continue;
index 87162a1..18df675 100644 (file)
@@ -25,14 +25,9 @@ ClassImp(AliAnalysisTaskScale)
 
 //________________________________________________________________________
 AliAnalysisTaskScale::AliAnalysisTaskScale() : 
-  AliAnalysisTaskSE(), 
-  fTracksName(), 
-  fClustersName(), 
-  fMinTrackPt(0.15),
-  fMinClusterPt(0.15),
+  AliAnalysisTaskEmcal("AliAnalysisTaskScale", kTRUE), 
   fScaleFunction(0),
   fGeom(0),
-  fOutputList(0), 
   fHistCentrality(0), 
   fHistPtTPCvsCent(0), 
   fHistPtEMCALvsCent(0), 
@@ -54,14 +49,9 @@ AliAnalysisTaskScale::AliAnalysisTaskScale() :
 
 //________________________________________________________________________
 AliAnalysisTaskScale::AliAnalysisTaskScale(const char *name) :
-  AliAnalysisTaskSE(name), 
-  fTracksName("Tracks"),
-  fClustersName("CaloClusters"),
-  fMinTrackPt(0.15),
-  fMinClusterPt(0.15),
+  AliAnalysisTaskEmcal(name, kTRUE), 
   fScaleFunction(0),
   fGeom(0),
-  fOutputList(0), 
   fHistCentrality(0), 
   fHistPtTPCvsCent(0), 
   fHistPtEMCALvsCent(0), 
@@ -80,7 +70,6 @@ AliAnalysisTaskScale::AliAnalysisTaskScale(const char *name) :
 {
   // Constructor.
 
-  DefineOutput(1, TList::Class());
 }
 
 //________________________________________________________________________
@@ -89,8 +78,8 @@ void AliAnalysisTaskScale::UserCreateOutputObjects()
   // Create my user objects.
 
   OpenFile(1);
-  fOutputList = new TList();
-  fOutputList->SetOwner();
+  fOutput = new TList();
+  fOutput->SetOwner();
 
   fHistCentrality         = new TH1F("Centrality","Centrality",              101, -1, 100);
   fHistPtTPCvsCent        = new TH2F("PtTPCvsCent","rho vs cent",            101, -1, 100, 500,   0, 1000);
@@ -108,23 +97,23 @@ void AliAnalysisTaskScale::UserCreateOutputObjects()
   fHistTrackEtaPhi        = new TH2F("TrackEtaPhi","Track eta phi",          100, -1.0, 1.0, 64,  0, 6.4);
   fHistClusterEtaPhi      = new TH2F("ClusterEtaPhi","Cluster eta phi",      100, -1.0, 1.0, 64, -3.2, 3.2);
 
-  fOutputList->Add(fHistCentrality);
-  fOutputList->Add(fHistPtTPCvsCent);
-  fOutputList->Add(fHistPtEMCALvsCent);
-  fOutputList->Add(fHistEtvsCent);
-  fOutputList->Add(fHistScalevsCent);
-  fOutputList->Add(fHistDeltaScalevsCent);
-  fOutputList->Add(fHistPtTPCvsNtrack);
-  fOutputList->Add(fHistPtEMCALvsNtrack);
-  fOutputList->Add(fHistEtvsNtrack);
-  fOutputList->Add(fHistScalevsNtrack);
-  fOutputList->Add(fHistDeltaScalevsNtrack);
-  fOutputList->Add(fHistTrackPtvsCent);
-  fOutputList->Add(fHistClusterPtvsCent);
-  fOutputList->Add(fHistTrackEtaPhi);
-  fOutputList->Add(fHistClusterEtaPhi);
-
-  PostData(1, fOutputList);
+  fOutput->Add(fHistCentrality);
+  fOutput->Add(fHistPtTPCvsCent);
+  fOutput->Add(fHistPtEMCALvsCent);
+  fOutput->Add(fHistEtvsCent);
+  fOutput->Add(fHistScalevsCent);
+  fOutput->Add(fHistDeltaScalevsCent);
+  fOutput->Add(fHistPtTPCvsNtrack);
+  fOutput->Add(fHistPtEMCALvsNtrack);
+  fOutput->Add(fHistEtvsNtrack);
+  fOutput->Add(fHistScalevsNtrack);
+  fOutput->Add(fHistDeltaScalevsNtrack);
+  fOutput->Add(fHistTrackPtvsCent);
+  fOutput->Add(fHistClusterPtvsCent);
+  fOutput->Add(fHistTrackEtaPhi);
+  fOutput->Add(fHistClusterEtaPhi);
+
+  PostData(1, fOutput);
 }
 
 //________________________________________________________________________
@@ -139,48 +128,25 @@ Double_t AliAnalysisTaskScale::GetScaleFactor(Double_t cent)
 }
 
 //________________________________________________________________________
-void AliAnalysisTaskScale::UserExec(Option_t *) 
+void AliAnalysisTaskScale::Init() 
 {
-  // Execute on each event.
+  // Init the analysis.
+
+  AliAnalysisTaskEmcal::Init();
 
-  if (!fGeom)
-    fGeom = AliEMCALGeometry::GetInstance();
+  fGeom = AliEMCALGeometry::GetInstance();
 
   if (!fGeom) {
     AliFatal("Can not create geometry");
     return;
   }
+}
 
-  // get input collections
-  TClonesArray *tracks = 0;
-  TClonesArray *clusters = 0;
-  TList *l = InputEvent()->GetList();
-  tracks = dynamic_cast<TClonesArray*>(l->FindObject(fTracksName));
-  if (!tracks) {
-    AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() ));
-    return;
-  }
-  clusters = dynamic_cast<TClonesArray*>(l->FindObject(fClustersName));
-  if (!clusters){
-    AliError(Form("Pointer to clusters %s == 0", fClustersName.Data() ));
-    return;
-  }
-
-  // get centrality
-  Double_t cent = -1;
-  AliCentrality *centrality = InputEvent()->GetCentrality();
-  if (centrality)
-    cent = centrality->GetCentralityPercentile("V0M");
-  if (cent < 0) {
-    AliError(Form("Centrality negative: %f", cent));
-    return;
-  }
-  fHistCentrality->Fill(cent);
+//________________________________________________________________________
+Bool_t AliAnalysisTaskScale::FillHistograms() 
+{
+  // Execute on each event.
 
-  // get vertex
-  Double_t vertex[3] = {0, 0, 0};
-  InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
-  
   const Double_t EmcalMinEta = fGeom->GetArm1EtaMin();
   const Double_t EmcalMaxEta = fGeom->GetArm1EtaMax();
   const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
@@ -195,68 +161,67 @@ void AliAnalysisTaskScale::UserExec(Option_t *)
   Double_t ptTPC   = 0;
   Double_t ptEMCAL = 0;
 
-  const Int_t Ntracks = tracks->GetEntries();
+  const Int_t Ntracks = fTracks->GetEntries();
   for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) {
-    AliVTrack *track = static_cast<AliVTrack*>(tracks->At(iTracks));
+    AliVTrack *track = static_cast<AliVTrack*>(fTracks->At(iTracks));
 
     if (!track)
       continue;
 
+    if (!AcceptTrack(track))
+      continue;
+
     if (TMath::Abs(track->Eta()) > 0.7)   // only accept tracks in the EMCal eta range
       continue;
 
-    fHistTrackPtvsCent->Fill(cent,track->Pt());
+    fHistTrackPtvsCent->Fill(fCent,track->Pt());
     fHistTrackEtaPhi->Fill(track->Eta(),track->Phi());
 
-    if (track->Pt()< fMinTrackPt)
-      continue;
-
     ptTPC += track->Pt();
     if ((track->Phi() > EmcalMaxPhi) || (track->Phi() < EmcalMinPhi))
       continue;
 
     ptEMCAL += track->Pt();
   }
+
+  if (ptTPC == 0) 
+    return kFALSE;
   
   Double_t Et = 0;
-  const Int_t Nclus = clusters->GetEntries();
+  const Int_t Nclus = fCaloClusters->GetEntries();
   for (Int_t iClus = 0; iClus < Nclus; ++iClus) {
-    AliVCluster *c = static_cast<AliVCluster*>(clusters->At(iClus));
+    AliVCluster *c = static_cast<AliVCluster*>(fCaloClusters->At(iClus));
     if (!c)
       continue;
-    if (!c->IsEMCAL())
+
+    if (!AcceptCluster(c))
       continue;
-    TLorentzVector nPart;
-    c->GetMomentum(nPart, vertex);
 
-    fHistClusterPtvsCent->Fill(cent,nPart.Pt());
-    fHistClusterEtaPhi->Fill(nPart.Eta(),nPart.Phi());
+    TLorentzVector nPart;
+    c->GetMomentum(nPart, fVertex);
 
-    if (nPart.Pt()< fMinClusterPt)
-      continue;
+    fHistClusterPtvsCent->Fill(fCent, nPart.Pt());
+    fHistClusterEtaPhi->Fill(nPart.Eta(), nPart.Phi());
 
     Et += nPart.Pt();
   }
-  
-  if (ptTPC == 0) 
-    return;
 
   const Double_t scalecalc = ((Et + ptEMCAL) / EmcalArea) * (TpcArea / ptTPC);
-
-  const Double_t scale     = GetScaleFactor(cent);
-
-  fHistPtTPCvsCent->Fill(cent, ptTPC);
-  fHistPtEMCALvsCent->Fill(cent, ptEMCAL);
-  fHistEtvsCent->Fill(cent, Et);
-  fHistScalevsCent->Fill(cent, scalecalc);
-  fHistDeltaScalevsCent->Fill(cent, scalecalc - scale);
+  const Double_t scale     = GetScaleFactor(fCent);
+
+  fHistCentrality->Fill(fCent);
+  fHistPtTPCvsCent->Fill(fCent, ptTPC);
+  fHistPtEMCALvsCent->Fill(fCent, ptEMCAL);
+  fHistEtvsCent->Fill(fCent, Et);
+  fHistScalevsCent->Fill(fCent, scalecalc);
+  fHistDeltaScalevsCent->Fill(fCent, scalecalc - scale);
   fHistPtTPCvsNtrack->Fill(Ntracks, ptTPC);
   fHistPtEMCALvsNtrack->Fill(Ntracks, ptEMCAL);
   fHistEtvsNtrack->Fill(Ntracks, Et);
   fHistScalevsNtrack->Fill(Ntracks, scalecalc);
   fHistDeltaScalevsNtrack->Fill(Ntracks, scalecalc - scale);
 
-  PostData(1, fOutputList);
+  return kTRUE;
 }      
 
 //________________________________________________________________________
index b53ec1e..8024f91 100644 (file)
@@ -9,35 +9,28 @@ class TH2F;
 class TF1;
 class AliEMCALGeometry;
 
-#include "AliAnalysisTaskSE.h"
+#include "AliAnalysisTaskEmcal.h"
 
-class AliAnalysisTaskScale : public AliAnalysisTaskSE {
+class AliAnalysisTaskScale : public AliAnalysisTaskEmcal {
  public:
   AliAnalysisTaskScale();
   AliAnalysisTaskScale(const char *name);
   virtual ~AliAnalysisTaskScale() {}
   
   virtual void           UserCreateOutputObjects();
-  virtual void           UserExec(Option_t *option);
   virtual void           Terminate(Option_t *);
 
-  void                   SetClustersName(const char *n)                        { fClustersName  = n    ; }
-  void                   SetMinClusterPt(Double_t min)                         { fMinClusterPt  = min  ; }
-  void                   SetMinTrackPt(Double_t min)                           { fMinTrackPt    = min  ; }
   void                   SetScaleFunction(TF1* sf)                             { fScaleFunction = sf   ; }
-  void                   SetTracksName(const char *n)                          { fTracksName    = n    ; }
   
  protected:
   virtual Double_t       GetScaleFactor(Double_t cent);
+  virtual Bool_t         FillHistograms();
+  void                   Init();
 
  private:
-  TString                fTracksName;             // name of track collection
-  TString                fClustersName;           // name of clusters collection
-  Double_t               fMinTrackPt;             // pt cut for scale factor calculation
-  Double_t               fMinClusterPt;           // pt cut for scale factor calculation
   TF1                   *fScaleFunction;          // scale factor as a function of centrality
+
   AliEMCALGeometry      *fGeom;                   //!ptr to emcal geometry object
-  TList                 *fOutputList;             //!output list
   TH1F                  *fHistCentrality;         //!output histogram
   TH2F                  *fHistPtTPCvsCent;        //!output histogram
   TH2F                  *fHistPtEMCALvsCent;      //!output histogram
@@ -57,6 +50,6 @@ class AliAnalysisTaskScale : public AliAnalysisTaskSE {
   AliAnalysisTaskScale(const AliAnalysisTaskScale&); // not implemented
   AliAnalysisTaskScale& operator=(const AliAnalysisTaskScale&); // not implemented
   
-  ClassDef(AliAnalysisTaskScale, 6); // Scale task
+  ClassDef(AliAnalysisTaskScale, 7); // Scale task
 };
 #endif
index 6c16bc2..63aac29 100644 (file)
@@ -71,6 +71,7 @@ class AliEmcalJet : public AliVParticle
   AliEmcalJet*      MatchedJet()                 const { return fMatched < 2 ? fClosestJets[fMatched] : 0; }
   Double_t          MaxClusterPt()               const { return MaxNeutralPt();            }
   Double_t          MaxTrackPt()                 const { return MaxChargedPt();            }
+  Double_t          MaxPartPt()                  const { return fMaxCPt < fMaxNPt ? fMaxNPt : fMaxCPt;     }
 
   void              AddClusterAt(Int_t clus, Int_t idx){ fClusterIDs.AddAt(clus, idx);     }
   void              AddTrackAt(Int_t track, Int_t idx) { fTrackIDs.AddAt(track, idx);      }
index e241529..448b4a6 100644 (file)
@@ -234,7 +234,6 @@ void AliEmcalJetTask::FindJets(TObjArray *tracks, TObjArray *clus, Int_t algo, D
       fjw.AddInputVector(nPart.Px(), nPart.Py(), nPart.Pz(), nPart.P(), -iClus - 100);  
     }
   }
-
   // run jet finder
   fjw.Run();
 
index c0f7640..c27b8e6 100644 (file)
@@ -49,12 +49,16 @@ void AliJetEmbeddingTask::Run()
   // Embed particles.
   
   if (fNClusters > 0 && fOutClusters) {
+    if (fCopyArray) 
+      CopyClusters();
     for (Int_t i = 0; i < fNClusters; ++i) {
       AddCluster();
     }
   }
  
   if (fNTracks > 0 && fOutTracks) {
+    if (fCopyArray) 
+      CopyTracks();
     for (Int_t i = 0; i < fNTracks; ++i) {
       AddTrack();
     }
index 8a10920..275c979 100644 (file)
@@ -7,10 +7,13 @@
 #include <TClonesArray.h>
 #include <TLorentzVector.h>
 #include <TRandom3.h>
+#include <TF1.h>
 
 #include "AliAnalysisManager.h"
 #include "AliVEvent.h"
 #include "AliVCluster.h"
+#include "AliESDCaloCluster.h"
+#include "AliAODCaloCluster.h"
 #include "AliEMCALDigit.h"
 #include "AliEMCALRecPoint.h"
 #include "AliPicoTrack.h"
@@ -39,6 +42,7 @@ AliJetModelBaseTask::AliJetModelBaseTask() :
   fCopyArray(kTRUE),
   fNClusters(0),
   fNTracks(0),
+  fPtSpectrum(0),
   fGeom(0),
   fClusters(0),
   fOutClusters(0),
@@ -66,6 +70,7 @@ AliJetModelBaseTask::AliJetModelBaseTask(const char *name) :
   fCopyArray(kTRUE),
   fNClusters(0),
   fNTracks(1),
+  fPtSpectrum(0),
   fGeom(0),
   fClusters(0),
   fOutClusters(0),
@@ -101,15 +106,16 @@ void AliJetModelBaseTask::Init()
     fPhiMax = fPhiMin;
   }
 
-  if (fNTracks > 0) {
+  if (fNTracks > 0 && !fTracks) {
     fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName));
     if (!fTracks) {
-      AliError(Form("Couldn't retrieve tracks with name %s!", fTracksName.Data()));
+      AliError(Form("%s: Couldn't retrieve tracks with name %s!", GetName(), fTracksName.Data()));
       return;
     }
-
-    if (strcmp(fTracks->GetClass()->GetName(), "AliPicoTrack")) {
-      AliError("Can only embed PicoTracks!");
+    
+    if (!fTracks->GetClass()->GetBaseClass("AliPicoTrack")) {
+      AliError(Form("%s: Collection %s does not contain AliPicoTrack objects!", GetName(), fTracksName.Data())); 
+      fTracks = 0;
       return;
     }
 
@@ -117,7 +123,7 @@ void AliJetModelBaseTask::Init()
       fOutTracksName = fTracksName;
       if (fCopyArray) {
        fOutTracksName += fSuffix;
-       fOutTracks = new TClonesArray(*fTracks);
+       fOutTracks = new TClonesArray("AliPicoTrack", fTracks->GetSize());
        fOutTracks->SetName(fOutTracksName);
       }
       else {
@@ -128,15 +134,21 @@ void AliJetModelBaseTask::Init()
     if (fCopyArray) {
       if (!(InputEvent()->FindListObject(fOutTracksName)))
        InputEvent()->AddObject(fOutTracks);
-      fOutTracks->Clear();
+      //fOutTracks->Clear();
     }
   }
 
-  if (fNClusters > 0) {
+  if (fNClusters > 0 && !fClusters) {
     fClusters = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName));
  
     if (!fClusters) {
-      AliError(Form("Couldn't retrieve clusters with name %s!", fCaloName.Data()));
+      AliError(Form("%s: Couldn't retrieve clusters with name %s!", GetName(), fCaloName.Data()));
+      return;
+    }
+
+    if (!fClusters->GetClass()->GetBaseClass("AliVCluster")) {
+      AliError(Form("%s: Collection %s does not contain AliVCluster objects!", GetName(), fCaloName.Data())); 
+      fClusters = 0;
       return;
     }
 
@@ -144,7 +156,7 @@ void AliJetModelBaseTask::Init()
       fOutCaloName = fCaloName;
       if (fCopyArray) {
        fOutCaloName += fSuffix;
-       fOutClusters = new TClonesArray(*fClusters);
+       fOutClusters = new TClonesArray(fClusters->GetClass()->GetName(), fClusters->GetSize());
        fOutClusters->SetName(fOutCaloName);
       }
       else {
@@ -155,7 +167,7 @@ void AliJetModelBaseTask::Init()
     if (fCopyArray) {
       if (!(InputEvent()->FindListObject(fOutCaloName)))
        InputEvent()->AddObject(fOutClusters);
-      fOutClusters->Clear();
+      //fOutClusters->Clear();
     }
 
     if (!fGeom) {
@@ -185,6 +197,13 @@ void AliJetModelBaseTask::Init()
     if (fPhiMin > EmcalMaxPhi) fPhiMin = EmcalMaxPhi;
     if (fPhiMin < EmcalMinPhi) fPhiMin = EmcalMinPhi;
   }
+
+  if (fCopyArray) {
+    if (fOutTracks)
+      fOutTracks->Clear();
+    if (fOutClusters)
+      fOutClusters->Clear();
+  }
 }
 
 //________________________________________________________________________
@@ -235,7 +254,10 @@ Double_t AliJetModelBaseTask::GetRandomPt()
 {
   // Get random pt.
 
-  return gRandom->Rndm() * (fPtMax - fPtMin) + fPtMin;
+  if (fPtSpectrum)
+    return fPtSpectrum->GetRandom();
+  else
+    return gRandom->Rndm() * (fPtMax - fPtMin) + fPtMin;
 }
 
 //________________________________________________________________________
@@ -323,13 +345,13 @@ AliPicoTrack* AliJetModelBaseTask::AddTrack(Double_t pt, Double_t eta, Double_t
     phi = GetRandomPhi();
 
   AliPicoTrack *track = new ((*fOutTracks)[nTracks]) AliPicoTrack(pt, 
-                                               eta, 
-                                               phi, 
-                                               1, 
-                                               100,    // MC flag!      
-                                               0, 
-                                               0, 
-                                               kFALSE);
+                                                                 eta, 
+                                                                 phi, 
+                                                                 1, 
+                                                                 100,    // MC flag!      
+                                                                 0, 
+                                                                 0, 
+                                                                 kFALSE);
   return track;
 }
 
@@ -340,6 +362,48 @@ void AliJetModelBaseTask::Run()
 }
 
 //________________________________________________________________________
+void AliJetModelBaseTask::CopyClusters()
+{
+  // Copy all the clusters in the new collection
+
+  Bool_t esdMode = (Bool_t)(fClusters->GetClass()->GetBaseClass("AliESDCaloCluster") != 0);
+  const Int_t nClusters = fClusters->GetEntriesFast();
+  
+  if (esdMode) {
+    for (Int_t i = 0; i < nClusters; ++i) {
+      AliESDCaloCluster *esdcluster = static_cast<AliESDCaloCluster*>(fClusters->At(i));
+      if (!esdcluster)
+       continue;
+      if (!esdcluster->IsEMCAL())
+       continue;
+      new ((*fOutClusters)[i]) AliESDCaloCluster(*esdcluster);
+    }
+  }
+  else {
+    for (Int_t i = 0; i < nClusters; ++i) {
+      AliAODCaloCluster *aodcluster = static_cast<AliAODCaloCluster*>(fClusters->At(i));
+      if (!aodcluster)
+       continue;
+      if (!aodcluster->IsEMCAL())
+       continue;
+      new ((*fOutClusters)[i]) AliAODCaloCluster(*aodcluster);
+    }
+  }
+}
+
+//________________________________________________________________________
+void AliJetModelBaseTask::CopyTracks()
+{
+  const Int_t nTracks = fTracks->GetEntriesFast();
+  for (Int_t i = 0; i < nTracks; ++i) {
+    AliPicoTrack *track = static_cast<AliPicoTrack*>(fTracks->At(i));
+    if (!track)
+      continue;
+    new ((*fOutTracks)[i]) AliPicoTrack(*track);
+  }
+}
+
+//________________________________________________________________________
 void AliJetModelBaseTask::UserExec(Option_t *) 
 {
   // Execute per event.
index 379171a..8f29b39 100644 (file)
@@ -7,6 +7,7 @@ class TClonesArray;
 class AliEMCALGeometry;
 class AliVCluster;
 class AliPicoTrack;
+class TF1;
 
 #include "AliAnalysisTaskSE.h"
 
@@ -30,6 +31,7 @@ class AliJetModelBaseTask : public AliAnalysisTaskSE {
   void         SetNTracks(Int_t n)                   { fNTracks      = n;    }
   void         SetGeometryName(const char *n)        { fGeomName     = n;    }
   void         SetSuffix(const char *s)              { fSuffix       = s;    }
+  void         SetPtSpectrum(TF1 *f)                 { fPtSpectrum   = f;    }
 
  protected:
 
@@ -41,6 +43,8 @@ class AliJetModelBaseTask : public AliAnalysisTaskSE {
   AliVCluster           *AddCluster(Double_t e = -1, Double_t eta = -999, Double_t phi = -1);   // add a cluster; if values are -1 generate random parameters
   AliVCluster           *AddCluster(Double_t e, Int_t absId);                                   // add a cluster with given energy and position
   AliPicoTrack          *AddTrack(Double_t pt = -1, Double_t eta = -999, Double_t phi = -1);    // add a track; if values are -1 generate random parameters
+  void                   CopyClusters();
+  void                   CopyTracks();
 
   TString                fGeomName;               // EMCal geometry name
   TString                fTracksName;             // name of track collection
@@ -57,6 +61,7 @@ class AliJetModelBaseTask : public AliAnalysisTaskSE {
   Bool_t                 fCopyArray;              // whether or not the array will be copied to a new one before modelling
   Int_t                  fNClusters;              // how many clusters are being processed
   Int_t                  fNTracks;                // how many tracks are being processed
+  TF1                   *fPtSpectrum;             // pt spectrum parametrization to extract random pt values
   AliEMCALGeometry      *fGeom;                   //!pointer to EMCal geometry
   TClonesArray          *fClusters;               //!cluster collection
   TClonesArray          *fOutClusters;            //!output cluster collection
@@ -67,6 +72,6 @@ class AliJetModelBaseTask : public AliAnalysisTaskSE {
   AliJetModelBaseTask(const AliJetModelBaseTask&);            // not implemented
   AliJetModelBaseTask &operator=(const AliJetModelBaseTask&); // not implemented
 
-  ClassDef(AliJetModelBaseTask, 1) // Jet modelling task
+  ClassDef(AliJetModelBaseTask, 2) // Jet modelling task
 };
 #endif
index bb36226..6607b83 100644 (file)
@@ -25,7 +25,7 @@ ClassImp(AliJetResponseMaker)
 
 //________________________________________________________________________
 AliJetResponseMaker::AliJetResponseMaker() : 
-  AliAnalysisTaskEmcalJet("AliJetResponseMaker"),
+  AliAnalysisTaskEmcalJet("AliJetResponseMaker", kTRUE),
   fMCTracksName("MCParticles"),
   fMCJetsName("MCJets"),
   fMaxDistance(0.25),
@@ -53,7 +53,7 @@ AliJetResponseMaker::AliJetResponseMaker() :
 
 //________________________________________________________________________
 AliJetResponseMaker::AliJetResponseMaker(const char *name) : 
-  AliAnalysisTaskEmcalJet(name),
+  AliAnalysisTaskEmcalJet(name, kTRUE),
   fMCTracksName("MCParticles"),
   fMCJetsName("MCJets"),
   fMaxDistance(0.25),
@@ -194,7 +194,7 @@ Bool_t AliJetResponseMaker::FillHistograms()
 
   for (Int_t i = 0; i < nMCJets; i++) {
 
-    AliEmcalJet* jet = dynamic_cast<AliEmcalJet*>(fMCJets->At(i));
+    AliEmcalJet* jet = static_cast<AliEmcalJet*>(fMCJets->At(i));
 
     if (!jet) {
       AliError(Form("Could not receive jet %d", i));
@@ -248,7 +248,7 @@ Bool_t AliJetResponseMaker::FillHistograms()
 
   for (Int_t i = 0; i < nJets; i++) {
 
-    AliEmcalJet* jet = dynamic_cast<AliEmcalJet*>(fJets->At(i));
+    AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(i));
 
     if (!jet) {
       AliError(Form("Could not receive mc jet %d", i));
@@ -298,7 +298,7 @@ void AliJetResponseMaker::DoJetLoop(TClonesArray *jets1, TClonesArray *jets2, Bo
 
   for (Int_t i = 0; i < nJets1; i++) {
 
-    AliEmcalJet* jet1 = dynamic_cast<AliEmcalJet*>(jets1->At(i));
+    AliEmcalJet* jet1 = static_cast<AliEmcalJet*>(jets1->At(i));
 
     if (!jet1) {
       AliError(Form("Could not receive jet %d", i));
@@ -307,14 +307,10 @@ void AliJetResponseMaker::DoJetLoop(TClonesArray *jets1, TClonesArray *jets2, Bo
 
     if (!AcceptJet(jet1, kTRUE, mc))
       continue;
-    
-    if (jet1->MaxTrackPt() < fPtBiasJetTrack && 
-        (fAnaType == kTPC || jet1->IsMC() || jet1->MaxClusterPt() < fPtBiasJetClus))
-       continue;
 
     for (Int_t j = 0; j < nJets2; j++) {
       
-      AliEmcalJet* jet2 = dynamic_cast<AliEmcalJet*>(jets2->At(j));
+      AliEmcalJet* jet2 = static_cast<AliEmcalJet*>(jets2->At(j));
       
       if (!jet2) {
        AliError(Form("Could not receive jet %d", j));
@@ -324,10 +320,6 @@ void AliJetResponseMaker::DoJetLoop(TClonesArray *jets1, TClonesArray *jets2, Bo
       if (!AcceptJet(jet2, kTRUE, !mc))
        continue;
       
-      if (jet2->MaxTrackPt() < fPtBiasJetTrack && 
-          (fAnaType == kTPC || jet2->IsMC() || jet2->MaxClusterPt() < fPtBiasJetClus))
-       continue;
-      
       Double_t deta = jet2->Eta() - jet1->Eta();
       Double_t dphi = jet2->Phi() - jet1->Phi();
       Double_t d = TMath::Sqrt(deta * deta + dphi * dphi);
@@ -351,17 +343,32 @@ Bool_t AliJetResponseMaker::RetrieveEventObjects()
   if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
     return kFALSE;
   
-  if (!fMCJetsName.IsNull()) {
+  if (!fMCJetsName.IsNull() && !fMCJets) {
     fMCJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fMCJetsName));
     if (!fMCJets) {
-      AliWarning(Form("Could not retrieve MC jets %s!", fMCJetsName.Data())); 
+      AliError(Form("%s: Could not retrieve mc jets %s!", GetName(), fMCJetsName.Data()));
+      return kFALSE;
+    }
+    else if (!fMCJets->GetClass()->GetBaseClass("AliEmcalJet")) {
+      AliError(Form("%s: Collection %s does not contain AliEmcalJet objects!", GetName(), fMCJetsName.Data())); 
+      fMCJets = 0;
+      return kFALSE;
     }
   }
 
-  if (!fMCTracksName.IsNull()) {
+  if (!fMCTracksName.IsNull() && !fMCTracks) {
     fMCTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fMCTracksName));
-    if (!fMCJets) {
-      AliWarning(Form("Could not retrieve MC tracks %s!", fMCTracksName.Data())); 
+    if (!fMCTracks) {
+      AliError(Form("%s: Could not retrieve mc tracks %s!", GetName(), fMCTracksName.Data())); 
+      return kFALSE;
+    }
+    else {
+      TClass *cl = fMCTracks->GetClass();
+      if (!cl->GetBaseClass("AliVParticle") && !cl->GetBaseClass("AliEmcalParticle")) {
+       AliError(Form("%s: Collection %s does not contain AliVParticle nor AliEmcalParticle objects!", GetName(), fMCTracksName.Data())); 
+       fMCTracks = 0;
+       return kFALSE;
+      }
     }
   }
 
index 0c6e5bd..3827338 100644 (file)
@@ -30,7 +30,8 @@ AliAnalysisTaskScale* AddTaskScale(
   TString name(Form("Scale_%s", nClusters));
   AliAnalysisTaskScale *scaletask = new AliAnalysisTaskScale(name);
   scaletask->SetTracksName(nTracks);
-  scaletask->SetClustersName(nClusters);
+  scaletask->SetClusName(nClusters);
+  scaletask->SetAnaType(AliAnalysisTaskEmcal::kEMCAL);
 
   //-------------------------------------------------------
   // Final settings, pass to manager and set the containers
index 51be78f..0561a72 100644 (file)
@@ -36,7 +36,8 @@ AliAnalysisTaskEmcal::AliAnalysisTaskEmcal() :
   fNbins(500),
   fMinBinPt(0),
   fMaxBinPt(250),
-  fPtCut(0.15),
+  fClusPtCut(0.15),
+  fTrackPtCut(0.15),
   fTracks(0),
   fCaloClusters(0),
   fCent(0),
@@ -62,7 +63,8 @@ AliAnalysisTaskEmcal::AliAnalysisTaskEmcal(const char *name, Bool_t histo) :
   fNbins(500),
   fMinBinPt(0),
   fMaxBinPt(250),
-  fPtCut(0.15),
+  fClusPtCut(0.15),
+  fTrackPtCut(0.15),
   fTracks(0),
   fCaloClusters(0),
   fCent(0),
@@ -92,9 +94,12 @@ void AliAnalysisTaskEmcal::UserExec(Option_t *)
 {
   // Main loop, called for each event.
 
-  if (!fInitialized) 
+  if (!fInitialized)
     Init();
 
+  if (!fInitialized)
+    return;
+
   if (!RetrieveEventObjects())
     return;
 
@@ -127,7 +132,7 @@ Bool_t AliAnalysisTaskEmcal::AcceptCluster(AliVCluster *clus, Bool_t acceptMC) c
   TLorentzVector nPart;
   clus->GetMomentum(nPart, const_cast<Double_t*>(fVertex));
 
-  if (nPart.Et() < fPtCut)
+  if (nPart.Et() < fClusPtCut)
     return kFALSE;
 
   return kTRUE;
@@ -144,7 +149,7 @@ Bool_t AliAnalysisTaskEmcal::AcceptEmcalPart(AliEmcalParticle *part, Bool_t acce
   if (fAnaType == kEMCAL && !part->IsEMCAL())
     return kFALSE;
 
-  if (part->Pt() < fPtCut)
+  if ((part->IsTrack() && part->Pt() < fTrackPtCut) || (part->IsCluster() && part->Pt() < fClusPtCut))
     return kFALSE;
 
   if (!acceptMC && part->IsMC())
@@ -164,7 +169,7 @@ Bool_t AliAnalysisTaskEmcal::AcceptTrack(AliVTrack *track, Bool_t acceptMC) cons
   if (!acceptMC && track->GetLabel() == 100)
     return kFALSE;
 
-  if (track->Pt() < fPtCut)
+  if (track->Pt() < fTrackPtCut)
     return kFALSE;
   
   return kTRUE;
@@ -205,20 +210,6 @@ void AliAnalysisTaskEmcal::Init()
 {
   // Init the analysis.
 
-  if (!fCaloName.IsNull() && (fAnaType == kEMCAL) && !fCaloClusters) {
-    fCaloClusters =  dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName));
-    if (!fCaloClusters) {
-      AliWarning(Form("%s: Could not retrieve clusters %s!", GetName(), fCaloName.Data())); 
-    }
-  }
-
-  if (!fTracksName.IsNull() && !fTracks) {
-    fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName));
-    if (!fTracks) {
-      AliWarning(Form("%s: Could not retrieve tracks %s!", GetName(), fTracksName.Data())); 
-    }
-  }
-
   SetInitialized();
 }
 
@@ -232,6 +223,38 @@ Bool_t AliAnalysisTaskEmcal::RetrieveEventObjects()
     return kFALSE;
   }
 
+  if (!fCaloName.IsNull() && (fAnaType == kEMCAL) && !fCaloClusters) {
+    fCaloClusters =  dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName));
+    if (!fCaloClusters) {
+      AliError(Form("%s: Could not retrieve clusters %s!", GetName(), fCaloName.Data())); 
+      return kFALSE;
+    }
+    else {
+      TClass *cl = fCaloClusters->GetClass();
+      if (!cl->GetBaseClass("AliVCluster") && !cl->GetBaseClass("AliEmcalParticle")) {
+       AliError(Form("%s: Collection %s does not contain AliVCluster nor AliEmcalParticle objects!", GetName(), fCaloName.Data())); 
+       fCaloClusters = 0;
+       return kFALSE;
+      }
+    }
+  }
+
+  if (!fTracksName.IsNull() && !fTracks) {
+    fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName));
+    if (!fTracks) {
+      AliError(Form("%s: Could not retrieve tracks %s!", GetName(), fTracksName.Data())); 
+      return kFALSE;
+    }
+    else {
+      TClass *cl = fTracks->GetClass();
+      if (!cl->GetBaseClass("AliVParticle") && !cl->GetBaseClass("AliEmcalParticle")) {
+       AliError(Form("%s: Collection %s does not contain AliVParticle nor AliEmcalParticle objects!", GetName(), fTracksName.Data())); 
+       fTracks = 0;
+       return kFALSE;
+      }
+    }
+  }
+
   fVertex[0] = 0;
   fVertex[1] = 0;
   fVertex[2] = 0;
index 9228812..cf3bb11 100644 (file)
@@ -37,7 +37,9 @@ class AliAnalysisTaskEmcal : public AliAnalysisTaskSE {
   void                        SetAnaType(EmcalAnaType type)                        { fAnaType        = type       ; }
   void                        SetClusName(const char *n)                           { fCaloName       = n          ; }
   void                        SetHistoBins(Int_t nbins, Float_t min, Float_t max)  { fNbins = nbins; fMinBinPt = min; fMaxBinPt = max; }
-  void                        SetPtCut(Float_t cut)                                { fPtCut          = cut        ; }
+  void                        SetPtCut(Float_t cut)                                { SetClusPtCut(cut); SetTrackPtCut(cut); }
+  void                        SetClusPtCut(Float_t cut)                            { fClusPtCut      = cut        ; }
+  void                        SetTrackPtCut(Float_t cut)                           { fTrackPtCut     = cut        ; }
   void                        SetTracksName(const char *n)                         { fTracksName     = n          ; }
 
  protected:
@@ -45,11 +47,11 @@ class AliAnalysisTaskEmcal : public AliAnalysisTaskSE {
   Bool_t                      AcceptEmcalPart(AliEmcalParticle *part,  Bool_t acceptMC = kFALSE)       const;
   Bool_t                      AcceptTrack(AliVTrack            *track, Bool_t acceptMC = kFALSE)       const;
   BeamType                    GetBeamType()                                                                 ;
-  virtual void                Init();
+  void                        Init();
   virtual Bool_t              FillHistograms()                                     { return fCreateHisto; }
-  void                        SetInitialized(Bool_t ini = kTRUE)                   { fInitialized    = ini        ; }
   virtual Bool_t              RetrieveEventObjects();
   virtual Bool_t              Run()                                                { return kTRUE                 ; }
+  void                        SetInitialized(Bool_t ini = kTRUE)                   { fInitialized    = ini        ; }
 
   EmcalAnaType                fAnaType;                    // analysis type
   Bool_t                      fInitialized;                // whether or not the task has been already initialized
@@ -59,7 +61,8 @@ class AliAnalysisTaskEmcal : public AliAnalysisTaskSE {
   Int_t                       fNbins;                      // no. of pt bins
   Float_t                     fMinBinPt;                   // min pt in histograms
   Float_t                     fMaxBinPt;                   // max pt in histograms
-  Float_t                     fPtCut;                      // cut on particle pt
+  Float_t                     fClusPtCut;                  // cut on cluster pt
+  Float_t                     fTrackPtCut;                 // cut on track pt
   TClonesArray               *fTracks;                     //!tracks
   TClonesArray               *fCaloClusters;               //!clusters
   Float_t                     fCent;                       //!event centrality
index abd3562..b712ff3 100644 (file)
@@ -77,10 +77,16 @@ void AliEmcalPicoTrackMaker::UserExec(Option_t *)
   }
 
   // retrieve tracks from input.
-  fTracksIn = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksInName));
-  if (!fTracksIn) {
-    AliError(Form("Could not retrieve tracks %s!", fTracksInName.Data())); 
-    return;
+  if (!fTracksIn) { 
+    fTracksIn = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksInName));
+    if (!fTracksIn) {
+      AliError(Form("Could not retrieve tracks %s!", fTracksInName.Data())); 
+      return;
+    }
+    if (!fTracksIn->GetClass()->GetBaseClass("AliVParticle")) {
+      AliError(Form("%s: Collection %s does not contain AliVParticle objects!", GetName(), fTracksInName.Data())); 
+      return;
+    }
   }
 
   // add tracks to event if not yet there
@@ -100,7 +106,7 @@ void AliEmcalPicoTrackMaker::UserExec(Option_t *)
   const Int_t Ntracks = fTracksIn->GetEntriesFast();
   for (Int_t iTracks = 0, nacc = 0; iTracks < Ntracks; ++iTracks) {
 
-    AliVTrack *track = dynamic_cast<AliVTrack*>(fTracksIn->At(iTracks));
+    AliVTrack *track = static_cast<AliVTrack*>(fTracksIn->At(iTracks));
 
     if (!track)
       continue;