]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWGJE/EMCALJetTasks/AliAnalysisTaskEmcalJetSample.cxx
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / AliAnalysisTaskEmcalJetSample.cxx
index 7a4226b6d23c67b4f8c3d39f724c5a622fef3f13..503bd06f1ee09ec4f97b05ba7a70b436b43dca4e 100644 (file)
@@ -7,15 +7,21 @@
 #include <TClonesArray.h>
 #include <TH1F.h>
 #include <TH2F.h>
+#include <TH3F.h>
 #include <TList.h>
 #include <TLorentzVector.h>
 
 #include "AliVCluster.h"
+#include "AliAODCaloCluster.h"
+#include "AliESDCaloCluster.h"
 #include "AliVTrack.h"
 #include "AliEmcalJet.h"
 #include "AliRhoParameter.h"
 #include "AliLog.h"
 #include "AliJetContainer.h"
+#include "AliParticleContainer.h"
+#include "AliClusterContainer.h"
+#include "AliPicoTrack.h"
 
 #include "AliAnalysisTaskEmcalJetSample.h"
 
@@ -23,12 +29,34 @@ ClassImp(AliAnalysisTaskEmcalJetSample)
 
 //________________________________________________________________________
 AliAnalysisTaskEmcalJetSample::AliAnalysisTaskEmcalJetSample() : 
-  AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetSample", kTRUE)
+  AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetSample", kTRUE),
+  fHistTracksPt(0),
+  fHistClustersPt(0),
+  fHistLeadingJetPt(0),
+  fHistJetsPhiEta(0),
+  fHistJetsPtArea(0),
+  fHistJetsPtLeadHad(0),
+  fHistJetsCorrPtArea(0),
+  fHistPtDEtaDPhiTrackClus(0),
+  fHistPtDEtaDPhiClusTrack(0),
+  fHistClustDx(0),
+  fHistClustDz(0),
+  fJetsCont(0),
+  fTracksCont(0),
+  fCaloClustersCont(0)
 
 {
   // Default constructor.
 
-  for (Int_t i = 0; i < 4; i++) {
+  fHistTracksPt       = new TH1*[fNcentBins];
+  fHistClustersPt     = new TH1*[fNcentBins];
+  fHistLeadingJetPt   = new TH1*[fNcentBins];
+  fHistJetsPhiEta     = new TH2*[fNcentBins];
+  fHistJetsPtArea     = new TH2*[fNcentBins];
+  fHistJetsPtLeadHad  = new TH2*[fNcentBins];
+  fHistJetsCorrPtArea = new TH2*[fNcentBins];
+
+  for (Int_t i = 0; i < fNcentBins; i++) {
     fHistTracksPt[i] = 0;
     fHistClustersPt[i] = 0;
     fHistLeadingJetPt[i] = 0;
@@ -43,12 +71,33 @@ AliAnalysisTaskEmcalJetSample::AliAnalysisTaskEmcalJetSample() :
 
 //________________________________________________________________________
 AliAnalysisTaskEmcalJetSample::AliAnalysisTaskEmcalJetSample(const char *name) : 
-  AliAnalysisTaskEmcalJet(name, kTRUE)
+  AliAnalysisTaskEmcalJet(name, kTRUE),
+  fHistTracksPt(0),
+  fHistClustersPt(0),
+  fHistLeadingJetPt(0),
+  fHistJetsPhiEta(0),
+  fHistJetsPtArea(0),
+  fHistJetsPtLeadHad(0),
+  fHistJetsCorrPtArea(0),
+  fHistPtDEtaDPhiTrackClus(0),
+  fHistPtDEtaDPhiClusTrack(0),
+  fHistClustDx(0),
+  fHistClustDz(0),
+  fJetsCont(0),
+  fTracksCont(0),
+  fCaloClustersCont(0)
 {
   // Standard constructor.
 
+  fHistTracksPt       = new TH1*[fNcentBins];
+  fHistClustersPt     = new TH1*[fNcentBins];
+  fHistLeadingJetPt   = new TH1*[fNcentBins];
+  fHistJetsPhiEta     = new TH2*[fNcentBins];
+  fHistJetsPtArea     = new TH2*[fNcentBins];
+  fHistJetsPtLeadHad  = new TH2*[fNcentBins];
+  fHistJetsCorrPtArea = new TH2*[fNcentBins];
 
-  for (Int_t i = 0; i < 4; i++) {
+  for (Int_t i = 0; i < fNcentBins; i++) {
     fHistTracksPt[i] = 0;
     fHistClustersPt[i] = 0;
     fHistLeadingJetPt[i] = 0;
@@ -74,9 +123,20 @@ void AliAnalysisTaskEmcalJetSample::UserCreateOutputObjects()
 
   AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
 
+  fJetsCont           = GetJetContainer(0);
+  if(fJetsCont) { //get particles and clusters connected to jets
+    fTracksCont       = fJetsCont->GetParticleContainer();
+    fCaloClustersCont = fJetsCont->GetClusterContainer();
+  } else {        //no jets, just analysis tracks and clusters
+    fTracksCont       = GetParticleContainer(0);
+    fCaloClustersCont = GetClusterContainer(0);
+  }
+  if(fTracksCont) fTracksCont->SetClassName("AliVTrack");
+  if(fCaloClustersCont) fCaloClustersCont->SetClassName("AliVCluster");
+
   TString histname;
 
-  for (Int_t i = 0; i < 4; i++) {
+  for (Int_t i = 0; i < fNcentBins; i++) {
     if (fParticleCollArray.GetEntriesFast()>0) {
       histname = "fHistTracksPt_";
       histname += i;
@@ -135,6 +195,21 @@ void AliAnalysisTaskEmcalJetSample::UserCreateOutputObjects()
       }
     }
   }
+
+  histname = "fHistPtDEtaDPhiTrackClus";
+  fHistPtDEtaDPhiTrackClus = new TH3F(histname.Data(),Form("%s;#it{p}_{T}^{track};#Delta#eta;#Delta#varphi",histname.Data()),100,0.,100.,100,-0.1,0.1,100,-0.1,0.1);
+  fOutput->Add(fHistPtDEtaDPhiTrackClus);
+
+  histname = "fHistPtDEtaDPhiClusTrack";
+  fHistPtDEtaDPhiClusTrack = new TH3F(histname.Data(),Form("%s;#it{p}_{T}^{clus};#Delta#eta;#Delta#varphi",histname.Data()),100,0.,100.,100,-0.1,0.1,100,-0.1,0.1);
+  fOutput->Add(fHistPtDEtaDPhiClusTrack);
+
+  fHistClustDx = new TH1F("fHistClustDx","fHistClustDx;Dx",1000,0.,1.);
+  fOutput->Add(fHistClustDx);
+
+  fHistClustDz = new TH1F("fHistClustDz","fHistClustDz;Dz",1000,0.,1.);
+  fOutput->Add(fHistClustDz);
+
   PostData(1, fOutput); // Post data for ALL output slots > 0 here.
 }
 
@@ -143,78 +218,129 @@ Bool_t AliAnalysisTaskEmcalJetSample::FillHistograms()
 {
   // Fill histograms.
 
-  if (fTracks) {
-    const Int_t ntracks = fTracks->GetEntriesFast();
-    for (Int_t it = 0; it < ntracks; it++) {
-      AliVTrack *track = static_cast<AliVTrack*>(fTracks->At(it));
-
-      if (!track) {
-       AliError(Form("Could not receive track %d", it));
-       continue;
-      }
-     
-      if (!AcceptTrack(track))
-       continue;
-
+  if (fTracksCont) {
+    AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle(0)); 
+    while(track) {
       fHistTracksPt[fCentBin]->Fill(track->Pt()); 
+      track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
     }
   }
-  
-  if (fCaloClusters) {
-    const Int_t nclusters = fCaloClusters->GetEntriesFast();
-    
-    for (Int_t ic = 0; ic < nclusters; ic++) {
-      AliVCluster *cluster = static_cast<AliVCluster*>(fCaloClusters->At(ic));
-      
-      if (!cluster) {
-       AliError(Form("Could not receive cluster %d", ic));
-       continue;
-      }
 
+  if (fCaloClustersCont) {
+    AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0); 
+    while(cluster) {
       TLorentzVector nPart;
       cluster->GetMomentum(nPart, fVertex);
       fHistClustersPt[fCentBin]->Fill(nPart.Pt());
+      Double_t dx = cluster->GetTrackDx();
+      Double_t dz = cluster->GetTrackDz();
+      fHistClustDx->Fill(dx);
+      fHistClustDz->Fill(dz);
+      cluster = fCaloClustersCont->GetNextAcceptCluster();
     }
   }
 
-  if (fJets) {
-
-    fJets->Sort();
-
-    const Int_t njets = fJets->GetEntriesFast();
-    Bool_t leadJet = kFALSE;
-    for (Int_t ij = 0; ij < njets; ij++) {
-
-      AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(ij));
-      if (!jet) {
-       AliError(Form("Could not receive jet %d", ij));
-       continue;
-      }  
-
-      if (!AcceptJet(jet))
-       continue;
-
-      if (!leadJet) {
-       fHistLeadingJetPt[fCentBin]->Fill(jet->Pt());
-       leadJet = kTRUE;
-      }
+  if (fJetsCont) {
+    AliEmcalJet *jet = fJetsCont->GetNextAcceptJet(0); 
+    while(jet) {
 
       fHistJetsPtArea[fCentBin]->Fill(jet->Pt(), jet->Area());
       fHistJetsPhiEta[fCentBin]->Fill(jet->Eta(), jet->Phi());
 
-      Float_t ptLeading = GetLeadingHadronPt(jet);
+      Float_t ptLeading = fJetsCont->GetLeadingHadronPt(jet);
       fHistJetsPtLeadHad[fCentBin]->Fill(jet->Pt(), ptLeading);
 
-      if (fRho) {
-       Float_t corrPt = jet->Pt() - fRhoVal * jet->Area();
+      if (fHistJetsCorrPtArea[fCentBin]) {
+       Float_t corrPt = jet->Pt() - fJetsCont->GetRhoVal() * jet->Area();
        fHistJetsCorrPtArea[fCentBin]->Fill(corrPt, jet->Area());
       }
+      jet = fJetsCont->GetNextAcceptJet(); 
     }
+    
+    jet = fJetsCont->GetLeadingJet();
+    if(jet) fHistLeadingJetPt[fCentBin]->Fill(jet->Pt());
   }
 
+  CheckClusTrackMatching();
+
   return kTRUE;
 }
 
+//________________________________________________________________________
+void AliAnalysisTaskEmcalJetSample::CheckClusTrackMatching()
+{
+  
+  if(!fTracksCont || !fCaloClustersCont)
+    return;
+
+  Double_t deta = 999;
+  Double_t dphi = 999;
+
+  //Get closest cluster to track
+  AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle(0)); 
+  while(track) {
+    //Get matched cluster
+    Int_t emc1 = track->GetEMCALcluster();
+    if(fCaloClustersCont && emc1>=0) {
+      AliVCluster *clusMatch = fCaloClustersCont->GetCluster(emc1);
+      if(clusMatch) {
+       AliPicoTrack::GetEtaPhiDiff(track, clusMatch, dphi, deta);
+       fHistPtDEtaDPhiTrackClus->Fill(track->Pt(),deta,dphi);
+      }
+    }
+    track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
+  }
+  
+  //Get closest track to cluster
+  AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0); 
+  while(cluster) {
+    TLorentzVector nPart;
+    cluster->GetMomentum(nPart, fVertex);
+    fHistClustersPt[fCentBin]->Fill(nPart.Pt());
+    
+    //Get matched track
+    AliVTrack *mt = NULL;      
+    AliAODCaloCluster *acl = dynamic_cast<AliAODCaloCluster*>(cluster);
+    if(acl) {
+      if(acl->GetNTracksMatched()>1)
+       mt = static_cast<AliVTrack*>(acl->GetTrackMatched(0));
+    }
+    else {
+      AliESDCaloCluster *ecl = dynamic_cast<AliESDCaloCluster*>(cluster);
+      Int_t im = ecl->GetTrackMatchedIndex();
+      if(fTracksCont && im>=0) {
+       mt = static_cast<AliVTrack*>(fTracksCont->GetParticle(im));
+      }
+    }
+    if(mt) {
+      AliPicoTrack::GetEtaPhiDiff(mt, cluster, dphi, deta);
+      fHistPtDEtaDPhiClusTrack->Fill(nPart.Pt(),deta,dphi);
+      
+      /* //debugging
+        if(mt->IsEMCAL()) {
+        Int_t emc1 = mt->GetEMCALcluster();
+        Printf("current id: %d  emc1: %d",fCaloClustersCont->GetCurrentID(),emc1);
+        AliVCluster *clm = fCaloClustersCont->GetCluster(emc1);
+        AliPicoTrack::GetEtaPhiDiff(mt, clm, dphi, deta);
+        Printf("deta: %f dphi: %f",deta,dphi);
+        }
+      */
+    }
+    cluster = fCaloClustersCont->GetNextAcceptCluster();
+  }
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskEmcalJetSample::ExecOnce() {
+
+  AliAnalysisTaskEmcalJet::ExecOnce();
+
+  if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0;
+  if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0;
+  if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0;
+
+}
+
 //________________________________________________________________________
 Bool_t AliAnalysisTaskEmcalJetSample::Run()
 {