PWGGA/EMCALJetTasks/AliHadCorrTask.cxx

   1 // $Id$
   2
   3 #include "AliHadCorrTask.h"
   4 #include <TClonesArray.h>
   5 #include <TParticle.h>
   6 #include "AliEmcalJet.h"
   7 #include "AliAnalysisManager.h"
   8 #include "AliESDtrack.h"
   9 #include "AliFJWrapper.h"
  10 #include "AliESDCaloCluster.h"
  11 #include "TList.h"
  12 #include "TH1F.h"
  13 #include "TH2F.h"
  14 #include "TChain.h"
  15 #include <TLorentzVector.h>
  16 #include <AliCentrality.h>
  17 #include "AliPicoTrack.h"
  18
  19 ClassImp(AliHadCorrTask)
  20
  21 //________________________________________________________________________
  22 AliHadCorrTask::AliHadCorrTask(const char *name) :
  23   AliAnalysisTaskSE("AliHadCorrTask"),
  24   fTracksName("Tracks"),
  25   fCaloName("CaloClusters"),
  26   fMinPt(0.15),
  27   fOutputList(0),
  28   fHistEbefore(0),
  29   fHistEafter(0),
  30   fHistNclusvsCent(0),
  31   fHistNclusMatchvsCent(0),
  32   fOutCaloName("CaloClustersOut"),
  33   fOutClusters(0)
  34 {
  35
  36   for(Int_t i=0; i<4; i++){
  37     fHistMatchEtaPhi[i]=0x0;
  38     fHistMatchEvsP[i]=0x0;
  39     fHistMatchdRvsEP[i]=0x0;
  40   }
  41
  42   // Standard constructor.
  43   cout << "Constructor for HadCorrTask " << name <<endl;
  44   if (!name)
  45     return;
  46
  47   SetName(name);
  48   fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.";
  49
  50   DefineInput(0,TChain::Class());
  51   DefineOutput(1,TList::Class());
  52 }
  53
  54 AliHadCorrTask::AliHadCorrTask() :
  55   AliAnalysisTaskSE("AliHadCorrTask"),
  56   fTracksName("Tracks"),
  57   fCaloName("CaloClusters"),
  58   fOutClusters(0)
  59 {
  60   // Standard constructor.
  61
  62   fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.";
  63 }
  64
  65
  66 //________________________________________________________________________
  67 AliHadCorrTask::~AliHadCorrTask()
  68 {
  69   // Destructor
  70 }
  71
  72 //________________________________________________________________________
  73 void AliHadCorrTask::UserCreateOutputObjects()
  74 {
  75
  76   fOutClusters = new TClonesArray("AliESDCaloCluster");
  77   fOutClusters->SetName(fOutCaloName);
  78
  79   OpenFile(1);
  80   fOutputList = new TList();
  81   fOutputList->SetOwner();
  82
  83   char name[200];
  84
  85   for(int icent=0; icent<4; icent++){
  86     sprintf(name,"fHistMatchEtaPhi_%i",icent);
  87     fHistMatchEtaPhi[icent] = new TH2F(name,name,400,-0.2,0.2,1600,-0.8,0.8);
  88     sprintf(name,"fHistMatchEvsP_%i",icent);
  89     fHistMatchEvsP[icent]=new TH2F(name,name,400,0.,200.,1000,0.,10.);
  90     sprintf(name,"fHistMatchdRvsEP_%i",icent);
  91     fHistMatchdRvsEP[icent]=new TH2F(name,name,1000,0.,1.,1000,0.,10.);
  92
  93
  94   fOutputList->Add(fHistMatchEtaPhi[icent]);
  95   fOutputList->Add(fHistMatchEvsP[icent]);
  96   fOutputList->Add(fHistMatchdRvsEP[icent]);
  97
  98   }
  99   fHistNclusvsCent=new TH1F("Nclusvscent","NclusVsCent",100,0,100);
 100   fHistNclusMatchvsCent=new TH1F("NclusMatchvscent","NclusMatchVsCent",100,0,100);
 101   fHistEbefore=new TH1F("Ebefore","Ebefore",100,0,100);
 102   fHistEafter=new TH1F("Eafter","Eafter",100,0,100);
 103
 104   fOutputList->Add(fHistNclusMatchvsCent);
 105   fOutputList->Add(fHistNclusvsCent);
 106   fOutputList->Add(fHistEbefore);
 107   fOutputList->Add(fHistEafter);
 108
 109   PostData(1, fOutputList);
 110 }
 111
 112 //________________________________________________________________________
 113 void AliHadCorrTask::UserExec(Option_t *)
 114 {
 115
 116   if (!(InputEvent()->FindListObject(fOutCaloName)))
 117     InputEvent()->AddObject(fOutClusters);
 118   else fOutClusters->Delete();
 119
 120
 121   AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
 122   TClonesArray *tracks = 0;
 123   TClonesArray *clus   = 0;
 124   TList *l = InputEvent()->GetList();
 125
 126   AliCentrality *centrality = 0;
 127   centrality = dynamic_cast<AliCentrality*>(l->FindObject("Centrality"));
 128   Float_t fCent=-1;
 129   if(centrality)
 130     fCent = centrality->GetCentralityPercentile("V0M");
 131   else
 132     fCent=99;//probably pp data
 133
 134   if(fCent<0) return;
 135
 136   if (fTracksName == "Tracks")
 137     am->LoadBranch("Tracks");
 138   tracks = dynamic_cast<TClonesArray*>(l->FindObject(fTracksName));
 139   if (!tracks) {
 140     AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() ));
 141     return;
 142   }
 143   const Int_t Ntrks = tracks->GetEntries();
 144
 145   if (fCaloName == "CaloClusters")
 146     am->LoadBranch("CaloClusters");
 147   clus = dynamic_cast<TClonesArray*>(l->FindObject(fCaloName));
 148   if (!clus) {
 149     AliError(Form("Pointer to clus %s == 0", fCaloName.Data() ));
 150     return;
 151   }
 152
 153   Int_t centbin=-1;
 154   if(fCent>=0 && fCent<10) centbin=0;
 155   else if(fCent>=10 && fCent<30) centbin=1;
 156   else if(fCent>=30 && fCent<50) centbin=2;
 157   else if(fCent>=50 && fCent<100) centbin=3;
 158
 159   if (clus) {
 160     Double_t vertex[3] = {0, 0, 0};
 161     InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
 162     const Int_t Nclus = clus->GetEntries();
 163     for (Int_t iClus = 0, iN = 0, clusCount=0; iClus < Nclus; ++iClus) {
 164       AliVCluster *c = dynamic_cast<AliVCluster*>(clus->At(iClus));
 165       if (!c)
 166         continue;
 167       if (!c->IsEMCAL())
 168         continue;
 169       TLorentzVector nPart;
 170
 171       c->SetEmcCpvDistance(-1);
 172       c->SetTrackDistance(999,999);
 173       Double_t dEtaMin  = 1e9;
 174       Double_t dPhiMin  = 1e9;
 175       Int_t    imin     = -1;
 176       for(Int_t t = 0; t<Ntrks; ++t) {
 177         AliVTrack *track = dynamic_cast<AliVTrack*>(tracks->At(t));
 178         Double_t etadiff=999;
 179         Double_t phidiff=999;
 180         AliPicoTrack::GetEtaPhiDiff(track,c,phidiff,etadiff);
 181         Double_t dR = TMath::Sqrt(etadiff*etadiff+phidiff*phidiff);
 182         Double_t dRmin = TMath::Sqrt(dEtaMin*dEtaMin+dPhiMin*dPhiMin);
 183         if(dR > 25)
 184           continue;
 185         if(dR<dRmin){
 186           dEtaMin = etadiff;
 187           dPhiMin = phidiff;
 188           imin = t;
 189         }
 190       }
 191       c->SetEmcCpvDistance(imin);
 192       c->SetTrackDistance(dPhiMin, dEtaMin);
 193       c->GetMomentum(nPart, vertex);
 194       Double_t energy = nPart.P();
 195       if(energy<fMinPt) continue;
 196       if (imin>=0) {
 197         Double_t dPhiMin = c->GetTrackDx();
 198         Double_t dEtaMin = c->GetTrackDz();
 199         fHistMatchEtaPhi[centbin]->Fill(dEtaMin,dPhiMin);
 200       }
 201
 202       fHistNclusvsCent->Fill(fCent);
 203
 204       if (fHadCorr>0) {
 205         if (imin>=0) {
 206           Double_t dPhiMin = c->GetTrackDx();
 207           Double_t dEtaMin = c->GetTrackDz();
 208           Double_t dR=TMath::Sqrt(dEtaMin*dEtaMin+dPhiMin*dPhiMin);
 209
 210           AliVTrack *t = dynamic_cast<AliVTrack*>(tracks->At(imin));
 211           if (t) {
 212             if (t->Pt()<fMinPt)
 213               continue;
 214             if (t->P()>0) fHistMatchEvsP[centbin]->Fill(energy,energy/t->P());
 215             if (t->P()>0) fHistMatchdRvsEP[centbin]->Fill(dR,energy/t->P());
 216             fHistEbefore->Fill(fCent,energy);
 217             if (dPhiMin<0.05 && dEtaMin<0.025) { // pp cuts!!!
 218               energy -= fHadCorr*t->P();
 219               fHistNclusMatchvsCent->Fill(fCent);
 220             }
 221             if (energy<0)
 222               continue;
 223           }
 224         }
 225         fHistEafter->Fill(fCent,energy);
 226
 227       }//end had correction if
 228       if (energy>0){//Output new corrected clusters
 229         AliESDCaloCluster *oc = new ((*fOutClusters)[clusCount]) AliESDCaloCluster(*(dynamic_cast<AliESDCaloCluster*>(c)));
 230         oc->SetE(energy);
 231         clusCount++;
 232       }
 233     }
 234
 235
 236   }
 237 }
 238
 239
 240 //________________________________________________________________________
 241 void AliHadCorrTask::Terminate(Option_t *)
 242 {
 243
 244 }
 245
 246 //________________________________________________________________________