Updates from Salvatore
authorloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 21 May 2012 10:16:39 +0000 (10:16 +0000)
committerloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 21 May 2012 10:16:39 +0000 (10:16 +0000)
PWGGA/EMCALJetTasks/AliAnalysisTaskSAJF.cxx
PWGGA/EMCALJetTasks/AliAnalysisTaskSAJF.h

index 005a4e8..e0803d2 100644 (file)
@@ -11,7 +11,6 @@
 #include <TH2F.h>
 #include <TList.h>
 #include <TLorentzVector.h>
-#include <TParticle.h>
 #include <TRandom3.h>
 #include <TParameter.h>
 
@@ -39,7 +38,6 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() :
   fTracksName("Tracks"),
   fCaloName("CaloClusters"),
   fJetsName("Jets"),
-  fKtJetsName("KtJets"),
   fEmbJetsName("EmbJets"),
   fRhoName("Rho"),
   fNbins(500),
@@ -50,7 +48,6 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() :
   fTracks(0),
   fCaloClusters(0),
   fJets(0),
-  fKtJets(0),
   fEmbJets(0),
   fCent(0),
   fCentBin(-1),
@@ -74,10 +71,11 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() :
   for (Int_t i = 0; i < 4; i++) {
     fHistJetsPt[i] = 0;
     fHistCorrJetsPt[i] = 0;
-    fHistJetsPtCutPart[i] = 0;
-    fHistCorrJetsPtCutPart[i] = 0;
-    fHistJetsNEF[i] = 0;
-    fHistJetsZ[i] = 0;
+    fHistUnfoldedJetsPt[i] = 0;
+    fHistJetsPtNonBias[i] = 0;
+    fHistCorrJetsPtNonBias[i] = 0;
+    fHistJetsNEFvsPt[i] = 0;
+    fHistJetsZvsPt[i] = 0;
     fHistLeadingJetPt[i] = 0;
     fHistCorrLeadingJetPt[i] = 0;
     fHist2LeadingJetPt[i] = 0;
@@ -113,7 +111,6 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) :
   fTracksName("Tracks"),
   fCaloName("CaloClusters"),
   fJetsName("Jets"),
-  fKtJetsName("KtJets"),
   fEmbJetsName("EmbJets"),
   fRhoName("Rho"),
   fNbins(500),
@@ -124,7 +121,6 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) :
   fTracks(0),
   fCaloClusters(0),
   fJets(0),
-  fKtJets(0),
   fEmbJets(0),
   fCent(0),
   fCentBin(-1),
@@ -147,10 +143,11 @@ AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) :
   for (Int_t i = 0; i < 4; i++) {
     fHistJetsPt[i] = 0;
     fHistCorrJetsPt[i] = 0;
-    fHistJetsPtCutPart[i] = 0;
-    fHistCorrJetsPtCutPart[i] = 0;
-    fHistJetsNEF[i] = 0;
-    fHistJetsZ[i] = 0;
+    fHistUnfoldedJetsPt[i] = 0;
+    fHistJetsPtNonBias[i] = 0;
+    fHistCorrJetsPtNonBias[i] = 0;
+    fHistJetsNEFvsPt[i] = 0;
+    fHistJetsZvsPt[i] = 0;
     fHistLeadingJetPt[i] = 0;
     fHistCorrLeadingJetPt[i] = 0;
     fHist2LeadingJetPt[i] = 0;
@@ -184,6 +181,8 @@ AliAnalysisTaskSAJF::~AliAnalysisTaskSAJF()
 void AliAnalysisTaskSAJF::UserCreateOutputObjects()
 {
   // Create histograms
+
+  Float_t binWidth = (fMaxPt - fMinPt) / fNbins;
   
   OpenFile(1);
   fOutput = new TList();
@@ -250,44 +249,51 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
     histname = "fHistJetsPt_";
     histname += i;
     fHistJetsPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt, fMaxPt);
-    fHistJetsPt[i]->GetXaxis()->SetTitle("E [GeV]");
+    fHistJetsPt[i]->GetXaxis()->SetTitle("Pt [GeV/c]");
     fHistJetsPt[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistJetsPt[i]);
 
     histname = "fHistCorrJetsPt_";
     histname += i;
-    fHistCorrJetsPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt, fMaxPt);
-    fHistCorrJetsPt[i]->GetXaxis()->SetTitle("E [GeV]");
+    fHistCorrJetsPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxPt, fMaxPt);
+    fHistCorrJetsPt[i]->GetXaxis()->SetTitle("Pt [GeV/c]");
     fHistCorrJetsPt[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistCorrJetsPt[i]);
 
-    histname = "fHistJetsPtCutPart_";
+    histname = "fHistUnfoldedJetsPt_";
+    histname += i;
+    fHistUnfoldedJetsPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxPt, fMaxPt);
+    fHistUnfoldedJetsPt[i]->GetXaxis()->SetTitle("Pt [GeV/c]");
+    fHistUnfoldedJetsPt[i]->GetYaxis()->SetTitle("counts");
+    fOutput->Add(fHistUnfoldedJetsPt[i]);
+
+    histname = "fHistJetsPtNonBias_";
     histname += i;
-    fHistJetsPtCutPart[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt, fMaxPt);
-    fHistJetsPtCutPart[i]->GetXaxis()->SetTitle("E [GeV]");
-    fHistJetsPtCutPart[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistJetsPtCutPart[i]);
+    fHistJetsPtNonBias[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt, fMaxPt);
+    fHistJetsPtNonBias[i]->GetXaxis()->SetTitle("Pt [GeV/c]");
+    fHistJetsPtNonBias[i]->GetYaxis()->SetTitle("counts");
+    fOutput->Add(fHistJetsPtNonBias[i]);
 
-    histname = "fHistCorrJetsPtCutPart_";
+    histname = "fHistCorrJetsPtNonBias_";
     histname += i;
-    fHistCorrJetsPtCutPart[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt, fMaxPt);
-    fHistCorrJetsPtCutPart[i]->GetXaxis()->SetTitle("E [GeV]");
-    fHistCorrJetsPtCutPart[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistCorrJetsPtCutPart[i]);
+    fHistCorrJetsPtNonBias[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxPt, fMaxPt);
+    fHistCorrJetsPtNonBias[i]->GetXaxis()->SetTitle("Pt [GeV/c]");
+    fHistCorrJetsPtNonBias[i]->GetYaxis()->SetTitle("counts");
+    fOutput->Add(fHistCorrJetsPtNonBias[i]);
     
-    histname = "fHistJetsNEF_";
+    histname = "fHistJetsNEFvsPt_";
     histname += i;
-    fHistJetsNEF[i] = new TH1F(histname.Data(), histname.Data(), fNbins, 0, 1.2);
-    fHistJetsNEF[i]->GetXaxis()->SetTitle("NEF");
-    fHistJetsNEF[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistJetsNEF[i]);
+    fHistJetsNEFvsPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, 0, 1.2, fNbins * 1.5, -fMaxPt * 0.75, fMaxPt * 0.75);
+    fHistJetsNEFvsPt[i]->GetXaxis()->SetTitle("NEF");
+    fHistJetsNEFvsPt[i]->GetYaxis()->SetTitle("Pt [GeV/c]");
+    fOutput->Add(fHistJetsNEFvsPt[i]);
 
-    histname = "fHistJetsZ_";
+    histname = "fHistJetsZvsPt_";
     histname += i;
-    fHistJetsZ[i] = new TH1F(histname.Data(), histname.Data(), fNbins, 0, 1.2);
-    fHistJetsZ[i]->GetXaxis()->SetTitle("Z");
-    fHistJetsZ[i]->GetYaxis()->SetTitle("counts");
-    fOutput->Add(fHistJetsZ[i]);
+    fHistJetsZvsPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, 0, 1.2, fNbins * 1.5, -fMaxPt * 0.75, fMaxPt * 0.75);
+    fHistJetsZvsPt[i]->GetXaxis()->SetTitle("Z");
+    fHistJetsZvsPt[i]->GetYaxis()->SetTitle("Pt [GeV/c]");
+    fOutput->Add(fHistJetsZvsPt[i]);
 
     histname = "fHistLeadingJetPt_";
     histname += i;
@@ -305,7 +311,7 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
 
     histname = "fHistCorrLeadingJetPt_";
     histname += i;
-    fHistCorrLeadingJetPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt, fMaxPt);
+    fHistCorrLeadingJetPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxPt, fMaxPt);
     fHistCorrLeadingJetPt[i]->GetXaxis()->SetTitle("P_{T} [GeV]");
     fHistCorrLeadingJetPt[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistCorrLeadingJetPt[i]);
@@ -361,14 +367,14 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
 
     histname = "fHistDeltaPtRC_";
     histname += i;
-    fHistDeltaPtRC[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt - fMaxPt / 2, fMaxPt / 2);
+    fHistDeltaPtRC[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt - fMaxPt / 2 + binWidth / 2, fMinPt + fMaxPt / 2 + binWidth / 2);
     fHistDeltaPtRC[i]->GetXaxis()->SetTitle("#deltaP_{T} [GeV/c]");
     fHistDeltaPtRC[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistDeltaPtRC[i]);
 
     histname = "fHistDeltaPtRCExLJ_";
     histname += i;
-    fHistDeltaPtRCExLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt - fMaxPt / 2, fMaxPt / 2);
+    fHistDeltaPtRCExLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt - fMaxPt / 2 + binWidth / 2, fMinPt + fMaxPt / 2 + binWidth / 2);
     fHistDeltaPtRCExLJ[i]->GetXaxis()->SetTitle("#deltaP_{T} [GeV/c]");
     fHistDeltaPtRCExLJ[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistDeltaPtRCExLJ[i]);
@@ -403,7 +409,7 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
 
     histname = "fHistDeltaPtMedKtEmb_";
     histname += i;
-    fHistDeltaPtMedKtEmb[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt - fMaxPt / 2, fMaxPt / 2);
+    fHistDeltaPtMedKtEmb[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinPt - fMaxPt / 2 + binWidth / 2, fMinPt + fMaxPt / 2 + binWidth / 2);
     fHistDeltaPtMedKtEmb[i]->GetXaxis()->SetTitle("#deltaP_{T} [GeV/c]");
     fHistDeltaPtMedKtEmb[i]->GetYaxis()->SetTitle("counts");
     fOutput->Add(fHistDeltaPtMedKtEmb[i]);
@@ -415,7 +421,7 @@ void AliAnalysisTaskSAJF::UserCreateOutputObjects()
 //________________________________________________________________________
 void AliAnalysisTaskSAJF::RetrieveEventObjects()
 {
-  if (strcmp(fKtJetsName,"") && fAnaType == kEMCAL) {
+  if (strcmp(fCaloName,"") && fAnaType == kEMCAL) {
     fCaloClusters =  dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName));
     if (!fCaloClusters) {
       AliWarning(Form("Could not retrieve clusters!")); 
@@ -432,13 +438,6 @@ void AliAnalysisTaskSAJF::RetrieveEventObjects()
     AliWarning(Form("Could not retrieve jets!")); 
   }
 
-  if (strcmp(fKtJetsName,"")) {
-    fKtJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fKtJetsName));
-    if (!fKtJets) {
-      AliWarning(Form("Could not retrieve Kt jets!")); 
-    }
-  }
-
   if (strcmp(fEmbJetsName,"")) {
     fEmbJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fEmbJetsName));
     if (!fEmbJets) {
@@ -446,15 +445,15 @@ void AliAnalysisTaskSAJF::RetrieveEventObjects()
     }
   }
 
-  fCent = InputEvent()->GetCentrality();
-  if (fCent) {
-    Float_t cent = fCent->GetCentralityPercentile("V0M");
-    if      (cent >=  0 && cent <   10) fCentBin = 0;
-    else if (cent >= 10 && cent <   30) fCentBin = 1;
-    else if (cent >= 30 && cent <   50) fCentBin = 2;
-    else if (cent >= 50 && cent <= 100) fCentBin = 3; 
+  AliCentrality *aliCent = InputEvent()->GetCentrality();
+  if (aliCent) {
+    fCent = aliCent->GetCentralityPercentile("V0M");
+    if      (fCent >=  0 && fCent <   10) fCentBin = 0;
+    else if (fCent >= 10 && fCent <   30) fCentBin = 1;
+    else if (fCent >= 30 && fCent <   50) fCentBin = 2;
+    else if (fCent >= 50 && fCent <= 100) fCentBin = 3; 
     else {
-      AliWarning(Form("Negative centrality: %f. Assuming 99", cent));
+      AliWarning(Form("Negative centrality: %f. Assuming 99", fCent));
       fCentBin = 3;
     }
   }
@@ -463,13 +462,17 @@ void AliAnalysisTaskSAJF::RetrieveEventObjects()
     fCentBin = 3;
   }
 
-  TParameter<Double_t> *rho = dynamic_cast<TParameter<Double_t>*>(InputEvent()->FindListObject(fRhoName));
+  fRho = -1;
   
-  if (rho) {
-    fRho = rho->GetVal();
-  }
-  else {
-    AliWarning("No rho information found in the event. Will recalculate it by myself...");
+  if (strcmp(fRhoName,"")) {
+    TParameter<Double_t> *rho = dynamic_cast<TParameter<Double_t>*>(InputEvent()->FindListObject(fRhoName));
+  
+    if (rho) {
+      fRho = rho->GetVal();
+    }
+    else {
+      AliWarning("Could not retrieve rho information.");
+    }
   }
 }
 
@@ -528,24 +531,6 @@ Int_t AliAnalysisTaskSAJF::GetNumberOfJets() const
 }
 
 //________________________________________________________________________
-AliEmcalJet* AliAnalysisTaskSAJF::GetKtJet(const Int_t i) const
-{
-  if (fKtJets)
-    return dynamic_cast<AliEmcalJet*>(fKtJets->At(i));
-  else
-    return 0;
-}
-
-//________________________________________________________________________
-Int_t AliAnalysisTaskSAJF::GetNumberOfKtJets() const
-{
-  if (fKtJets)
-    return fKtJets->GetEntriesFast();
-  else
-    return 0;
-}
-
-//________________________________________________________________________
 AliEmcalJet* AliAnalysisTaskSAJF::GetEmbJet(const Int_t i) const
 {
   if (fEmbJets)
@@ -568,26 +553,6 @@ void AliAnalysisTaskSAJF::FillHistograms()
 {
   Float_t A = fJetRadius * fJetRadius * TMath::Pi();
 
-  Float_t cent = 100;
-
-  // as a rule will use rho provided by a previous RhoTask (see RetrieveEventObjects)
-  if (fRho < 0) {
-    Float_t rhoKt = 0;
-    DoKtJetLoop(rhoKt);
-    
-    if (rhoKt == 0) {
-      AliWarning("Skipping event with less than 2 reconstructed kt jets...");
-      return;
-    }
-
-    fRho = rhoKt;
-  }
-
-  if (fCent)
-    cent = fCent->GetCentralityPercentile("V0M");
-
-  fHistCentrality->Fill(cent);
-
   Int_t maxJetIndex  = -1;
   Int_t max2JetIndex = -1;
 
@@ -600,6 +565,8 @@ void AliAnalysisTaskSAJF::FillHistograms()
   if (!jet) 
     return;
 
+  fHistCentrality->Fill(fCent);
+
   fHistLeadingJetPt[fCentBin]->Fill(jet->Pt());
   jet->SortConstituents();
   
@@ -696,7 +663,7 @@ void AliAnalysisTaskSAJF::FillHistograms()
       maxEmbPartPhi = track->Phi();
     }
     else {
-      AliWarning("Embedded particle type not recognized (neither AliVCluster nor AliVTrack)!");
+      AliWarning(Form("%s - Embedded particle type not recognized (neither AliVCluster nor AliVTrack)!", GetName()));
       return;
     }
     fHistEmbJets[fCentBin]->Fill(maxJet->Pt());
@@ -708,7 +675,7 @@ void AliAnalysisTaskSAJF::FillHistograms()
   }
   else {
     if (maxPart != 0)
-      AliWarning("Embedded particle is not the leading particle of the leading jet!");
+      AliWarning(Form("%s - Embedded particle is not the leading particle of the leading jet!", GetName()));
   }
 }
 
@@ -739,22 +706,22 @@ void AliAnalysisTaskSAJF::DoJetLoop(Int_t &maxJetIndex, Int_t &max2JetIndex)
 
     Float_t corrPt = jet->Pt() - fRho * jet->Area();
 
-    fHistJetsPt[fCentBin]->Fill(jet->Pt());
+    fHistJetsPtNonBias[fCentBin]->Fill(jet->Pt());
+    fHistCorrJetsPtNonBias[fCentBin]->Fill(corrPt);
 
-    if (corrPt > 0)
-      fHistCorrJetsPt[fCentBin]->Fill(corrPt);
+    if (!AcceptJet(jet, kTRUE))
+      continue;
 
-    fHistJetPhiEta->Fill(jet->Eta(), jet->Phi());
-    fHistJetsNEF[fCentBin]->Fill(jet->NEF());
+    fHistJetsPt[fCentBin]->Fill(jet->Pt());
+    fHistCorrJetsPt[fCentBin]->Fill(corrPt);
 
-    Float_t maxPt = 0;
+    fHistJetPhiEta->Fill(jet->Eta(), jet->Phi());
+    fHistJetsNEFvsPt[fCentBin]->Fill(jet->NEF(), corrPt);
 
     for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
       AliVTrack *track = jet->TrackAt(it, fTracks);
       if (track)
-       fHistJetsZ[fCentBin]->Fill(track->Pt() / jet->Pt());
-      if (track->Pt() > maxPt)
-       maxPt = track->Pt();
+       fHistJetsZvsPt[fCentBin]->Fill(track->Pt() / jet->Pt(), corrPt);
     }
 
     for (Int_t ic = 0; ic < jet->GetNumberOfClusters(); ic++) {
@@ -763,114 +730,24 @@ void AliAnalysisTaskSAJF::DoJetLoop(Int_t &maxJetIndex, Int_t &max2JetIndex)
       if (cluster) {
        TLorentzVector nPart;
        cluster->GetMomentum(nPart, vertex);
-       fHistJetsZ[fCentBin]->Fill(nPart.Et() / jet->Pt());
-
-       if (nPart.Et() > maxPt)
-         maxPt = nPart.Et();
+       fHistJetsZvsPt[fCentBin]->Fill(nPart.Et() / jet->Pt(), corrPt);
       }
     }
 
-    if (maxPt >= fPtCutJetPart) {
-      fHistJetsPtCutPart[fCentBin]->Fill(jet->Pt());
-      if (corrPt > 0)
-       fHistCorrJetsPtCutPart[fCentBin]->Fill(corrPt);
-    }
-
-    if (maxJetPt < jet->Pt()) {
+    if (maxJetPt < corrPt) {
       max2JetPt = maxJetPt;
       max2JetIndex = maxJetIndex;
-      maxJetPt = jet->Pt();
+      maxJetPt = corrPt;
       maxJetIndex = ij;
     }
-    else if (max2JetPt < jet->Pt()) {
-      max2JetPt = jet->Pt();
+    else if (max2JetPt < corrPt) {
+      max2JetPt = corrPt;
       max2JetIndex = ij;
     }
   } //jet loop 
 }
 
 //________________________________________________________________________
-void AliAnalysisTaskSAJF::DoKtJetLoop(Float_t &rhoKt, Int_t nLJs)
-{
-  rhoKt = 0;
-
-  Int_t nktjets =  GetNumberOfKtJets();
-
-  Int_t NoOfZeroJets = 0;
-  if (nktjets > 0) {
-    
-    TArrayI ktJets(nktjets);
-    ktJets.Reset(-1);
-    for (Int_t ij = 0; ij < nktjets; ij++) {
-      
-      AliEmcalJet* jet = GetKtJet(ij);
-      
-      if (!jet) {
-       AliError(Form("Could not receive jet %d", ij));
-       continue;
-      } 
-      
-      if (jet->Pt() <= 0 || jet->Area() <= 0) {
-       NoOfZeroJets++;
-       continue;
-      }
-      
-      if (!AcceptJet(jet)) {
-       NoOfZeroJets++;
-       continue;
-      }
-    
-      Float_t rho = jet->Pt() / jet->Area();
-      Int_t i = nktjets;
-      AliEmcalJet *cmpjet = 0;
-      do {
-       i--;
-       if (ktJets[i] >= 0)
-         cmpjet = GetKtJet(ktJets[i]);
-       else
-         cmpjet = 0;
-      } while (cmpjet && rho < cmpjet->Pt() / cmpjet->Area() && i > 0);
-      
-      memmove(ktJets.GetArray() + nktjets - ij - 1, ktJets.GetArray() + nktjets - ij, (ij + i - nktjets + 1) * sizeof(Int_t));
-      ktJets[i] = ij;
-    } //kt jet loop 
-
-    nktjets -= NoOfZeroJets;
-    if (nktjets < 1) return;
-
-    memmove(ktJets.GetArray(), ktJets.GetArray() + NoOfZeroJets, nktjets * sizeof(Int_t));
-
-    nktjets -= nLJs;
-    if (nktjets < 1) return;
-
-    const Float_t maxEta = fMaxEta - fJetRadius;
-    const Float_t minEta = fMinEta + fJetRadius;
-    const Float_t maxPhi = fMaxPhi - fJetRadius;
-    const Float_t minPhi = fMinPhi + fJetRadius;
-
-    for (Int_t i = 0; i < nktjets; i++) {
-      AliEmcalJet *cmpjet = GetKtJet(ktJets[i]);
-      if (cmpjet->Eta() > maxEta || cmpjet->Eta() < minEta
-         || cmpjet->Phi() > maxPhi || cmpjet->Phi() < minPhi) {
-       nktjets--;
-       memmove(ktJets.GetArray() + i, ktJets.GetArray() + i + 1, (nktjets - i) * sizeof(Int_t));
-       i--;
-      }
-    }
-
-    if (nktjets % 2 != 0) {   // odd
-      AliEmcalJet *medJet = GetKtJet(ktJets[(nktjets - 1) / 2]);
-      rhoKt = medJet->Pt() / medJet->Area();
-    }
-    else {   // even
-      AliEmcalJet *medJet1 = GetKtJet(ktJets[nktjets / 2]);
-      AliEmcalJet *medJet2 = GetKtJet(ktJets[nktjets / 2 - 1]);
-      rhoKt = (medJet1->Pt() / medJet1->Area() + medJet2->Pt() / medJet2->Area()) / 2;
-    }
-  }
-}
-
-//________________________________________________________________________
 Bool_t AliAnalysisTaskSAJF::DoEmbJetLoop(AliEmcalJet* &maxJet, TObject* &maxPart)
 {
   Double_t vertex[3] = {0, 0, 0};
@@ -898,7 +775,7 @@ Bool_t AliAnalysisTaskSAJF::DoEmbJetLoop(AliEmcalJet* &maxJet, TObject* &maxPart
     if (jet->Pt() <= 0)
        continue;
  
-    if (!AcceptJet(jet))
+    if (!AcceptJet(jet, kTRUE))
       continue;
     
     if (jet->Pt() > maxJetPt) {
@@ -1196,8 +1073,11 @@ Bool_t AliAnalysisTaskSAJF::AcceptJet(Float_t eta, Float_t phi) const
 }
 
 //________________________________________________________________________
-Bool_t AliAnalysisTaskSAJF::AcceptJet(AliEmcalJet* jet) const
+Bool_t AliAnalysisTaskSAJF::AcceptJet(AliEmcalJet* jet, Bool_t checkPt) const
 {
+  if (checkPt && jet->MaxTrackPt() < fPtCutJetPart && jet->MaxClusterPt() < fPtCutJetPart)
+    return kFALSE;
+
   return AcceptJet(jet->Eta(), jet->Phi());
 }
 
index e385055..6ef5124 100644 (file)
@@ -31,16 +31,15 @@ class AliAnalysisTaskSAJF : public AliAnalysisTaskSE {
   void                        Terminate(Option_t *option);
   void                        Init();
 
-  void                        SetTracksName(const char *n)                         { fTracksName    = n          ; }
-  void                        SetClusName(const char *n)                           { fCaloName      = n          ; }
-  void                        SetJetsName(const char *n)                           { fJetsName      = n          ; }
-  void                        SetKtJetsName(const char *n)                         { fKtJetsName    = n          ; }
-  void                        SetEmbJetsName(const char *n)                        { fEmbJetsName   = n          ; }
-  void                        SetRhoName(const char *n)                            { fRhoName       = n          ; }
-  void                        SetAnaType(SAJFAnaType type)                         { fAnaType       = type       ; }
-  void                        SetJetRadius(Float_t r)                              { fJetRadius     = r          ; } 
-  void                        SetPtCut(Float_t cut)                                { fPtCut         = cut        ; }
-  void                        SetPtCutJetPart(Float_t cut)                         { fPtCutJetPart  = cut        ; }
+  void                        SetTracksName(const char *n)                         { fTracksName     = n          ; }
+  void                        SetClusName(const char *n)                           { fCaloName       = n          ; }
+  void                        SetJetsName(const char *n)                           { fJetsName       = n          ; }
+  void                        SetEmbJetsName(const char *n)                        { fEmbJetsName    = n          ; }
+  void                        SetRhoName(const char *n)                            { fRhoName        = n          ; }
+  void                        SetAnaType(SAJFAnaType type)                         { fAnaType        = type       ; }
+  void                        SetJetRadius(Float_t r)                              { fJetRadius      = r          ; } 
+  void                        SetPtCut(Float_t cut)                                { fPtCut          = cut        ; }
+  void                        SetPtCutJetPart(Float_t cut)                         { fPtCutJetPart   = cut        ; }
   void                        SetHistoBins(Int_t nbins, Float_t min, Float_t max)  { fNbins = nbins; fMinPt = min; fMaxPt = max; }
 
  protected:
@@ -51,13 +50,11 @@ class AliAnalysisTaskSAJF : public AliAnalysisTaskSE {
   Int_t                       GetNumberOfCaloClusters()                                            const;
   AliEmcalJet                *GetJet(const Int_t i)                                                const;
   Int_t                       GetNumberOfJets()                                                    const;
-  AliEmcalJet                *GetKtJet(const Int_t i)                                              const;
-  Int_t                       GetNumberOfKtJets()                                                  const;
   AliEmcalJet                *GetEmbJet(const Int_t i)                                             const;
   Int_t                       GetNumberOfEmbJets()                                                 const;
   Bool_t                      AcceptTrack(AliVTrack* track, Bool_t acceptMC = kFALSE)              const;
   Bool_t                      AcceptCluster(AliVCluster* clus, Bool_t acceptMC = kFALSE)                ;
-  Bool_t                      AcceptJet(AliEmcalJet* jet)                                          const;
+  Bool_t                      AcceptJet(AliEmcalJet* jet, Bool_t checkPt = kFALSE)                 const;
   Bool_t                      AcceptJet(Float_t eta, Float_t phi)                                  const;
   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;
@@ -66,13 +63,11 @@ class AliAnalysisTaskSAJF : public AliAnalysisTaskSE {
   void                        RetrieveEventObjects()                                                                             ;
   void                        FillHistograms()                                                                                   ;
   void                        DoJetLoop(Int_t &maxJetIndex, Int_t &max2JetIndex)                                                 ;
-  void                        DoKtJetLoop(Float_t &rhoKt, Int_t nLJs = 2)                                                        ;
   Bool_t                      DoEmbJetLoop(AliEmcalJet* &maxJet, TObject* &maxPart)                                              ;
   void                        DoTrackLoop(Float_t &rhoTracksLJex, Float_t &rhoTracks, Int_t maxJetIndex, Int_t max2JetIndex)     ;
   void                        DoClusterLoop(Float_t &rhoClusLJex, Float_t &rhoClus,Int_t maxJetIndex, Int_t max2JetIndex)        ;
   Bool_t                      GetRigidConePt(Float_t &pt, Float_t &eta, Float_t &phi, AliEmcalJet *jet = 0,  Float_t minD = 1.0) ;
 
-
   SAJFAnaType                 fAnaType;                    // Analysis type
   Float_t                     fMinEta;                     // Minimum eta accepatance
   Float_t                     fMaxEta;                     // Maximum eta accepatance
@@ -82,7 +77,6 @@ class AliAnalysisTaskSAJF : public AliAnalysisTaskSE {
   TString                     fTracksName;                 // Name of track collection
   TString                     fCaloName;                   // Name of calo cluster collection
   TString                     fJetsName;                   // Name of jet collection
-  TString                     fKtJetsName;                 // Name of kt jet collection
   TString                     fEmbJetsName;                // Name of embedded jets collection
   TString                     fRhoName;                    // Name of rho object
   Int_t                       fNbins;                      // No. of pt bins
@@ -93,9 +87,8 @@ class AliAnalysisTaskSAJF : public AliAnalysisTaskSE {
   TClonesArray               *fTracks;                     //!Tracks
   TClonesArray               *fCaloClusters;               //!Clusters
   TClonesArray               *fJets;                       //!Jets
-  TClonesArray               *fKtJets;                     //!Kt Jets
   TClonesArray               *fEmbJets;                    //!Embedded Jets
-  AliCentrality              *fCent;                       //!Event centrality
+  Float_t                     fCent;                       //!Event centrality
   Int_t                       fCentBin;                    //!Event centrality bin
   Float_t                     fRho;                        //!Event rho
   TList                      *fOutput;                     //!Output list
@@ -112,10 +105,11 @@ class AliAnalysisTaskSAJF : public AliAnalysisTaskSE {
   TH2F                       *fHistRCPhiEta;               //!Phi-Eta distribution of embedded jets
   TH1F                       *fHistJetsPt[4];              //!Inclusive jet pt spectrum
   TH1F                       *fHistCorrJetsPt[4];          //!Corrected inclusive jet pt spectrum
-  TH1F                       *fHistJetsPtCutPart[4];       //!Inclusive jet pt spectrum selecting jet with a particle minimum fPtCutJetPart
-  TH1F                       *fHistCorrJetsPtCutPart[4];   //!Corrected inclusive jet pt spectrum selecting jet with a particle minimum fPtCutJetPart
-  TH1F                       *fHistJetsNEF[4];             //!Jet neutral energy fraction
-  TH1F                       *fHistJetsZ[4];               //!Constituent Pt over Jet Pt ratio
+  TH1F                       *fHistUnfoldedJetsPt[4];      //!Unfolded inclusive jet pt spectrum
+  TH1F                       *fHistJetsPtNonBias[4];       //!Inclusive jet pt spectrum selecting jet with a particle minimum fPtCutJetPart
+  TH1F                       *fHistCorrJetsPtNonBias[4];   //!Corrected inclusive jet pt spectrum selecting jet with a particle minimum fPtCutJetPart
+  TH2F                       *fHistJetsNEFvsPt[4];         //!Jet neutral energy fraction
+  TH2F                       *fHistJetsZvsPt[4];           //!Constituent Pt over Jet Pt ratio
   TH1F                       *fHistLeadingJetPt[4];        //!Leading jet pt spectrum
   TH1F                       *fHistCorrLeadingJetPt[4];    //!Leading jet pt spectrum
   TH1F                       *fHist2LeadingJetPt[4];       //!Second leading jet pt spectrum