]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWGGA/EMCALJetTasks/AliAnalysisTaskScale.cxx
Changes for task:
[u/mrichter/AliRoot.git] / PWGGA / EMCALJetTasks / AliAnalysisTaskScale.cxx
index 1cbde06ffa6518e87a829b4a52e9f26ca0172533..fd037f0fbaa47aae572ed9d2a85c49bfaf3f423f 100644 (file)
 // $Id$
+//
+// Scale task.
+//
+// Author: R.Reed, M.Connors
+
+#include "AliAnalysisTaskScale.h"
 
-#include "TChain.h"
-#include "TTree.h"
-#include "TList.h"
-#include "TH1.h"
-#include "TH2.h"
 #include <TClonesArray.h>
-#include <TParticle.h>
+#include <TF1.h>
+#include <TH1F.h>
+#include <TH2F.h>
+#include <TList.h>
 #include <TLorentzVector.h>
 
-#include "AliAnalysisTask.h"
 #include "AliAnalysisManager.h"
-#include "AliVCluster.h"
-
-#include "AliESDEvent.h"
-#include "AliESDInputHandler.h"
 #include "AliCentrality.h"
-
-#include "AliAnalysisTaskScale.h"
+#include "AliEMCALGeometry.h"
+#include "AliLog.h"
+#include "AliVCluster.h"
+#include "AliVEvent.h"
+#include "AliVTrack.h"
 
 ClassImp(AliAnalysisTaskScale)
 
 //________________________________________________________________________
-AliAnalysisTaskScale::AliAnalysisTaskScale(const char *name) 
-  : AliAnalysisTaskSE(name), fESD(0), fOutputList(0), fHistCentrality(0), fHistPtTPCvsCent(0), fHistPtEMCALvsCent(0), fHistEtvsCent(0),  fHistScalevsCent(0),  fHistDeltaScalevsCent(0), fHistPtTPCvsNtrack(0), fHistPtEMCALvsNtrack(0), fHistEtvsNtrack(0),  fHistScalevsNtrack(0),  fHistDeltaScalevsNtrack(0),
-    fTracksName("tracks"),
-    fClustersName("clusters")
+AliAnalysisTaskScale::AliAnalysisTaskScale() : 
+  AliAnalysisTaskEmcal("AliAnalysisTaskScale", kTRUE), 
+  fScaleFunction(0),
+  fGeom(0),
+  fHistCentrality(0), 
+  fHistPtTPCvsCent(0), 
+  fHistPtEMCALvsCent(0), 
+  fHistEtvsCent(0),  
+  fHistScalevsCent(0),  
+  fHistDeltaScalevsCent(0), 
+  fHistPtTPCvsNtrack(0), 
+  fHistPtEMCALvsNtrack(0), 
+  fHistEtvsNtrack(0),  
+  fHistScalevsNtrack(0),  
+  fHistDeltaScalevsNtrack(0),
+  fHistTrackPtvsCent(0), 
+  fHistClusterPtvsCent(0),
+  fHistTrackEtaPhi(0), 
+  fHistClusterEtaPhi(0)
+{
+  // Default constructor.
+}
 
+//________________________________________________________________________
+AliAnalysisTaskScale::AliAnalysisTaskScale(const char *name) :
+  AliAnalysisTaskEmcal(name, kTRUE), 
+  fScaleFunction(0),
+  fGeom(0),
+  fHistCentrality(0), 
+  fHistPtTPCvsCent(0), 
+  fHistPtEMCALvsCent(0), 
+  fHistEtvsCent(0),  
+  fHistScalevsCent(0),  
+  fHistDeltaScalevsCent(0), 
+  fHistPtTPCvsNtrack(0), 
+  fHistPtEMCALvsNtrack(0), 
+  fHistEtvsNtrack(0),  
+  fHistScalevsNtrack(0),  
+  fHistDeltaScalevsNtrack(0),
+  fHistTrackPtvsCent(0), 
+  fHistClusterPtvsCent(0),
+  fHistTrackEtaPhi(0), 
+  fHistClusterEtaPhi(0)
 {
-  // Constructor
+  // Constructor.
 
-  DefineInput(0, TChain::Class());
-  DefineOutput(1, TList::Class());
 }
 
 //________________________________________________________________________
 void AliAnalysisTaskScale::UserCreateOutputObjects()
 {
-  OpenFile(1);
-  fOutputList = new TList();
-  fOutputList->SetOwner();
+  // Create my user objects.
 
-  fHistCentrality = new TH1F("fHistCentrality","Centrality",101,-1,100);
+  OpenFile(1);
+  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);
+  fHistPtEMCALvsCent      = new TH2F("PtEMCALvsCent","rho vs cent",          101, -1, 100, 500,   0, 1000);
+  fHistEtvsCent           = new TH2F("EtvsCent","rho vs cent",               101, -1, 100, 500,   0, 1000);
+  fHistScalevsCent        = new TH2F("ScalevsCent","rho vs cent",            101, -1, 100, 400,   0, 4);
+  fHistDeltaScalevsCent   = new TH2F("DeltaScalevsCent","rho vs cent",       101, -1, 100, 400,  -2, 2);
+  fHistPtTPCvsNtrack      = new TH2F("PtTPCvsNtrack","rho vs cent",          500,  0, 2500, 500,  0, 1000);
+  fHistPtEMCALvsNtrack    = new TH2F("PtEMCALvsNtrack","rho vs cent",        500,  0, 2500, 500,  0, 1000);
+  fHistEtvsNtrack         = new TH2F("EtvsNtrack","rho vs cent",             500,  0, 2500, 500,  0, 1000);
+  fHistScalevsNtrack      = new TH2F("ScalevsNtrack","rho vs cent",          500,  0, 2500, 400,  0, 4);
+  fHistDeltaScalevsNtrack = new TH2F("DeltaScalevsNtrack","rho vs cent",     500,  0, 2500, 400, -2, 2);
+  fHistTrackPtvsCent      = new TH2F("TrackPtvsCent","Track pt vs cent",     101, -1, 100,  500,  0, 100);
+  fHistClusterPtvsCent    = new TH2F("ClusterPtvsCent","Cluster pt vs cent", 101, -1, 100,  500,  0, 100);
+  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);
+
+  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);
+}
 
-  fHistPtTPCvsCent = new TH2F("fHistPtTPCvsCent","rho vs cent",101,-1,100,500,0,1000);
-  fHistPtEMCALvsCent = new TH2F("fHistPtEMCALvsCent","rho vs cent",101,-1,100,500,0,1000);
-  fHistEtvsCent = new TH2F("fHistEtvsCent","rho vs cent",101,-1,100,500,0,1000);
-  fHistScalevsCent = new TH2F("fHistScalevsCent","rho vs cent",101,-1,100,400,0,4.0);
-  fHistDeltaScalevsCent = new TH2F("fHistDeltaScalevsCent","rho vs cent",101,-1,100,400,-2.0,2.0);
-  
-  
-  fHistPtTPCvsNtrack = new TH2F("fHistPtTPCvsNtrack","rho vs cent",500,0,2500,500,0,1000);
-  fHistPtEMCALvsNtrack = new TH2F("fHistPtEMCALvsNtrack","rho vs cent",500,0,2500,500,0,1000);
-  fHistEtvsNtrack = new TH2F("fHistEtvsNtrack","rho vs cent",500,0,2500,500,0,1000);
-  fHistScalevsNtrack = new TH2F("fHistScalevsNtrack","rho vs cent",500,0,2500,400,0,4.0);
-  fHistDeltaScalevsNtrack = new TH2F("fHistDeltaScalevsNtrack","rho vs cent",500,0,2500,400,-2.0,2.0);
-    
-  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);
-  
+//________________________________________________________________________
+Double_t AliAnalysisTaskScale::GetScaleFactor(Double_t cent)
+{
+  // Get scale function.
 
-  PostData(1, fOutputList);
+  Double_t scale = -1;
+  if (fScaleFunction)
+    scale = fScaleFunction->Eval(cent);
+  return scale;
 }
 
 //________________________________________________________________________
-void AliAnalysisTaskScale::UserExec(Option_t *
+void AliAnalysisTaskScale::ExecOnce(
 {
-  fESD = dynamic_cast<AliESDEvent*>(InputEvent());
-  if (!fESD) {
-    printf("ERROR: fESD not available\n");
-    return;
-  }
+  // Init the analysis.
 
-  AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
-  TClonesArray *tracks = 0;
-  TClonesArray *clusters = 0;
-  AliCentrality *centrality = 0;
-  TList *l = InputEvent()->GetList();
+  fGeom = AliEMCALGeometry::GetInstance();
 
-  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() ));
+  if (!fGeom) {
+    AliFatal("Can not create geometry");
     return;
   }
-  
 
-  centrality = dynamic_cast<AliCentrality*>(l->FindObject("Centrality"));
-  float fCent = centrality->GetCentralityPercentile("V0M");
-  fHistCentrality->Fill(fCent);
+  AliAnalysisTaskEmcal::ExecOnce();
+}
+
+//________________________________________________________________________
+Bool_t AliAnalysisTaskScale::FillHistograms() 
+{
+  // Execute on each event.
+
+  const Double_t EmcalMinEta = fGeom->GetArm1EtaMin();
+  const Double_t EmcalMaxEta = fGeom->GetArm1EtaMax();
+  const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
+  const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
 
-  float scale = 0.000066*fCent*fCent-0.0015*fCent+1.5; //march 25th fit
-  float ptTPC = 0;float ptEMCAL=0;
-  const Int_t Ntracks = tracks->GetEntries();
+  const Double_t TpcMinPhi   = 0;
+  const Double_t TpcMaxPhi   = 2 * TMath::Pi();
+
+  const Double_t TpcArea     = (TpcMaxPhi - TpcMinPhi) * (EmcalMinEta - EmcalMaxEta);
+  const Double_t EmcalArea   = (EmcalMaxPhi - EmcalMinPhi) * (EmcalMinEta - EmcalMaxEta);
+
+  Double_t ptTPC   = 0;
+  Double_t ptEMCAL = 0;
+
+  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 (fabs(track->Eta())>0.7)
+
+    if (!AcceptTrack(track))
       continue;
-    ptTPC+=track->Pt();
-    if ((track->Phi()>3.14)||(track->Phi()<1.4))
+
+    if (TMath::Abs(track->Eta()) > 0.7)   // only accept tracks in the EMCal eta range
       continue;
-    ptEMCAL+=track->Pt();
-  } //track loop 
-  
-  Double_t vertex[3] = {0, 0, 0};
-  
-  float Et = 0;
-  InputEvent()->GetPrimaryVertex()->GetXYZ(vertex);
-  const Int_t Nclus = clusters->GetEntries();
-  for (Int_t iClus = 0, iN = 0; iClus < Nclus; ++iClus) {
-    
-    AliVCluster *c = dynamic_cast<AliVCluster*>(clusters->At(iClus));
-    if (!c->IsEMCAL())
+
+    fHistTrackPtvsCent->Fill(fCent,track->Pt());
+    fHistTrackEtaPhi->Fill(track->Eta(),track->Phi());
+
+    ptTPC += track->Pt();
+    if ((track->Phi() > EmcalMaxPhi) || (track->Phi() < EmcalMinPhi))
       continue;
-    TLorentzVector nPart;
-    c->GetMomentum(nPart, vertex);
-    Et += nPart.P();
+
+    ptEMCAL += track->Pt();
   }
+
+  if (ptTPC == 0) 
+    return kFALSE;
   
-  float scalecalc = ((Et+ptEMCAL)/2.44)*(8.8/ptTPC); //2.44 = EMCAL Acc, 8.8 = TPC Acc
-  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);
-}      
+  Double_t Et = 0;
+  const Int_t Nclus = fCaloClusters->GetEntries();
+  for (Int_t iClus = 0; iClus < Nclus; ++iClus) {
+    AliVCluster *c = static_cast<AliVCluster*>(fCaloClusters->At(iClus));
+    if (!c)
+      continue;
 
-//________________________________________________________________________
-void AliAnalysisTaskScale::Terminate(Option_t *) 
-{}
+    if (!AcceptCluster(c))
+      continue;
+
+    TLorentzVector nPart;
+    c->GetMomentum(nPart, fVertex);
+
+    fHistClusterPtvsCent->Fill(fCent, nPart.Pt());
+    fHistClusterEtaPhi->Fill(nPart.Eta(), nPart.Phi());
 
+    Et += nPart.Pt();
+  }
 
+  const Double_t scalecalc = ((Et + ptEMCAL) / EmcalArea) * (TpcArea / ptTPC);
+  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);
+
+  return kTRUE;
+}      
 
+//________________________________________________________________________
+void AliAnalysisTaskScale::Terminate(Option_t *) 
+{
+  // Called once at the end of the analysis.
+}