PWG4/totEt/AliAnalysisTaskHadEt.cxx

   1 //_________________________________________________________________________
   2 //  Utility Class for transverse energy studies; charged hadrons
   3 //  Task for analysis
   4 //  - reconstruction and MC output
   5 // implementation file
   6 //
   7 //Created by Christine Nattrass, Rebecca Scott, Irakli Martashvili
   8 //University of Tennessee at Knoxville
   9 //_________________________________________________________________________
  10 //Necessary to read config macros
  11 #include <TROOT.h>
  12 #include <TSystem.h>
  13 #include <TInterpreter.h>
  14
  15 #include "TChain.h"
  16 #include "TList.h"
  17 #include "TH2F.h"
  18
  19 #include "AliESDEvent.h"
  20 #include "AliMCEvent.h"
  21 #include "AliESDtrackCuts.h"
  22
  23 #include "AliAnalysisTaskHadEt.h"
  24 #include "AliAnalysisHadEtReconstructed.h"
  25 #include "AliAnalysisHadEtMonteCarlo.h"
  26
  27 #include <iostream>
  28
  29 using namespace std;
  30
  31 ClassImp(AliAnalysisTaskHadEt)
  32
  33
  34
  35 //________________________________________________________________________
  36   AliAnalysisTaskHadEt::AliAnalysisTaskHadEt(const char *name, Bool_t isMc, TString recoConfigFile, TString mcConfigFile) :
  37         AliAnalysisTaskTransverseEnergy(name, isMc)
  38         ,fRecAnalysis(0)
  39         ,fMCAnalysis(0)
  40         ,fIsSim(isMc)
  41 {
  42     // Constructor
  43   fMCConfigFile = mcConfigFile;
  44   fRecoConfigFile = recoConfigFile;
  45
  46   if(fMCAnalysis) delete fMCAnalysis;
  47   if(fRecAnalysis) delete fRecAnalysis;
  48
  49   if (fRecoConfigFile.Length()) {
  50     cout<<"Rereading AliAnalysisHadEtReconstructed configuration file "<<fRecoConfigFile<<endl;
  51     gROOT->LoadMacro(fRecoConfigFile);
  52     fRecAnalysis = (AliAnalysisHadEtReconstructed *) gInterpreter->ProcessLine("ConfigHadEtReconstructed()");
  53   }
  54
  55   if (fMCConfigFile.Length()) {
  56     cout<<"Rereading AliAnalysisHadEtMonteCarlo configuration file "<<fMCConfigFile<<endl;
  57     gROOT->LoadMacro(fMCConfigFile);
  58     fMCAnalysis = (AliAnalysisHadEtMonteCarlo *) gInterpreter->ProcessLine("ConfigHadEtMonteCarlo()");
  59     fMCAnalysis->SetHadEtReconstructed(fRecAnalysis);
  60   }
  61
  62     // Define input and output slots here
  63     // Input slot #0 works with a TChain
  64     DefineInput(0, TChain::Class());
  65     // Output slot #1 writes into a TH1 container
  66
  67     DefineOutput(1, TList::Class());
  68 }
  69 AliAnalysisTaskHadEt::~AliAnalysisTaskHadEt(){//Destructor
  70   delete fRecAnalysis;
  71   delete fMCAnalysis;
  72 }
  73
  74
  75 //________________________________________________________________________
  76 void AliAnalysisTaskHadEt::UserCreateOutputObjects()
  77 {
  78     // Create histograms
  79
  80     // Called once
  81
  82
  83
  84   fOutputList = new TList;
  85   fOutputList->SetOwner();
  86   fMCAnalysis->SetHistoList(fOutputList);
  87   fRecAnalysis->SetHistoList(fOutputList);
  88   if(fIsSim) fMCAnalysis->CreateHistograms();
  89   fRecAnalysis->CreateHistograms();
  90
  91
  92   if(fRecAnalysis->DataSet() != fMCAnalysis->DataSet()){
  93     cout<<"Warning: Reconstruction data set and Monte Carlo data set are not the same!  Setting data set to "<<fRecAnalysis->DataSet()<<endl;
  94   }
  95
  96   Bool_t selectPrimaries=kTRUE;
  97   if(fEsdtrackCutsITSTPC) delete fEsdtrackCutsITSTPC;
  98   if(fEsdtrackCutsITS) delete fEsdtrackCutsITS;
  99   if(fEsdtrackCutsTPC) delete fEsdtrackCutsTPC;
 100   //We do not use these because we are using the 2010 900 GeV data
 101 //   if(fRecAnalysis->DataSet()==2009){
 102 //     cout<<"Setting track cuts for the 2009 p+p collisions at 900 GeV"<<endl;
 103 //     fEsdtrackCutsITSTPC = AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(selectPrimaries);
 104 //     fEsdtrackCutsITSTPC->SetName("fEsdTrackCuts");
 105 //     fEsdtrackCutsTPC = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
 106 //     fEsdtrackCutsTPC->SetName("fEsdTrackCutsTPCOnly");
 107 //     //ITS stand alone cuts - similar to 2009 cuts but with only ITS hits required
 108 //     fEsdtrackCutsITS =  AliESDtrackCuts::GetStandardITSPureSATrackCuts2009(kTRUE,kFALSE);//we do want primaries but we do not want to require PID info
 109 //     fEsdtrackCutsITS->SetName("fEsdTrackCutsITS");
 110 //   }
 111   if(fRecAnalysis->DataSet()==2010 || fRecAnalysis->DataSet()==20111||fRecAnalysis->DataSet()==2009){
 112     if(fRecAnalysis->DataSet()==2010)cout<<"Setting track cuts for the 2010 p+p collisions at 7 TeV"<<endl;
 113     else{cout<<"Setting track cuts for the 2011 p+p collisions at 2.76 TeV"<<endl;}
 114     //cout<<"Warning:  Have not set 2010 track cuts yet!!"<<endl;
 115     fEsdtrackCutsITSTPC = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(selectPrimaries);
 116     fEsdtrackCutsITSTPC->SetName("fEsdTrackCuts");
 117     fEsdtrackCutsTPC = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
 118     fEsdtrackCutsTPC->SetName("fEsdTrackCutsTPCOnly");
 119     //ITS stand alone cuts - similar to 2009 cuts but with only ITS hits required
 120     fEsdtrackCutsITS =  AliESDtrackCuts::GetStandardITSPureSATrackCuts2010(kTRUE,kFALSE);//we do want primaries but we do not want to require PID info
 121     fEsdtrackCutsITS->SetName("fEsdTrackCutsITS");
 122   }
 123   if(fRecAnalysis->DataSet()==20100){
 124     cout<<"Setting track cuts for the 2010 Pb+Pb collisions at 2.76 TeV"<<endl;
 125     //cout<<"Warning:  Have not set 2010 track cuts yet!!"<<endl;
 126     fEsdtrackCutsITSTPC = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(selectPrimaries);
 127     fEsdtrackCutsITSTPC->SetName("fEsdTrackCuts");
 128     fEsdtrackCutsTPC = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
 129     fEsdtrackCutsTPC->SetName("fEsdTrackCutsTPCOnly");
 130     //ITS stand alone cuts - similar to 2009 cuts but with only ITS hits required
 131     fEsdtrackCutsITS =  AliESDtrackCuts::GetStandardITSSATrackCutsPbPb2010(kTRUE,kFALSE);//we do want primaries but we do not want to require PID info
 132     // fEsdtrackCutsITS =  AliESDtrackCuts::GetStandardITSPureSATrackCuts2010(kTRUE,kFALSE);//we do want primaries but we do not want to require PID info
 133    fEsdtrackCutsITS->SetName("fEsdTrackCutsITS");
 134   }
 135
 136   fOutputList->Add(fEsdtrackCutsITSTPC);
 137   fOutputList->Add(fEsdtrackCutsTPC);
 138   fOutputList->Add(fEsdtrackCutsITS);
 139   if(fEsdtrackCutsITSTPC && fEsdtrackCutsTPC){
 140     fRecAnalysis->SetITSTrackCuts( GetITSTrackCuts());
 141     fMCAnalysis->SetITSTrackCuts( GetITSTrackCuts());
 142     fRecAnalysis->SetTPCITSTrackCuts( GetTPCITSTrackCuts());
 143     fMCAnalysis->SetTPCITSTrackCuts( GetTPCITSTrackCuts());
 144     fRecAnalysis->SetTPCOnlyTrackCuts( GetTPCOnlyTrackCuts());
 145     fMCAnalysis->SetTPCOnlyTrackCuts( GetTPCOnlyTrackCuts());
 146     //add ITS stuff!
 147   }
 148   else{
 149     Printf("Error: no track cuts!");
 150   }
 151 }
 152
 153 //________________________________________________________________________
 154 void AliAnalysisTaskHadEt::UserExec(Option_t *)
 155 { // execute method
 156   fESDEvent = dynamic_cast<AliESDEvent*>(InputEvent());
 157 if (!fESDEvent) {
 158   Printf("ERROR: Could not retrieve event");
 159   return;
 160  }
 161 //cout<<"AliAnalysisHadEtReconstructed 90"<<endl;
 162
 163 Int_t res = CheckPhysicsSelection(fESDEvent->GetRunNumber()); // Check if the physics selection is valid for this run
 164
 165 AliCentrality *cent = GetCentralityObject();
 166
 167 if(res == 0 && cent){
 168
 169   //cout<<"New Event"<<endl;
 170
 171   fRecAnalysis->AnalyseEvent(fESDEvent);
 172
 173   AliMCEvent* mcEvent = MCEvent();
 174 // if (!mcEvent) {
 175 //   Printf("ERROR: Could not retrieve MC event");
 176 //  }
 177   if (mcEvent && fESDEvent && fIsSim){
 178       ((AliAnalysisHadEtMonteCarlo*)fMCAnalysis)->AnalyseEvent((AliVEvent*)mcEvent,(AliVEvent*)fESDEvent);
 179       if(fMCAnalysis->Full()){
 180         fMCAnalysis->FillSimTotEtMinusRecoTotEtFullAcceptanceTPC( fRecAnalysis->GetCorrectedTotEtFullAcceptanceTPC() );
 181         fMCAnalysis->FillSimTotEtMinusRecoTotEtFullAcceptanceITS( fRecAnalysis->GetCorrectedTotEtFullAcceptanceITS() );
 182         fMCAnalysis->FillSimTotEtMinusRecoTotEtFullAcceptanceTPCNoPID( fRecAnalysis->GetCorrectedTotEtFullAcceptanceTPCNoPID() );
 183         fMCAnalysis->FillSimTotEtMinusRecoTotEtFullAcceptanceITSNoPID( fRecAnalysis->GetCorrectedTotEtFullAcceptanceITSNoPID() );
 184         fMCAnalysis->FillSimHadEtMinusRecoHadEtFullAcceptanceTPC( fRecAnalysis->GetCorrectedHadEtFullAcceptanceTPC() );
 185         fMCAnalysis->FillSimHadEtMinusRecoHadEtFullAcceptanceITS( fRecAnalysis->GetCorrectedHadEtFullAcceptanceITS() );
 186         fMCAnalysis->FillSimHadEtMinusRecoHadEtFullAcceptanceTPCNoPID( fRecAnalysis->GetCorrectedHadEtFullAcceptanceTPCNoPID() );
 187         fMCAnalysis->FillSimHadEtMinusRecoHadEtFullAcceptanceITSNoPID( fRecAnalysis->GetCorrectedHadEtFullAcceptanceITSNoPID() );
 188       }
 189       if(fMCAnalysis->EMCAL()){
 190         fMCAnalysis->FillSimTotEtMinusRecoTotEtEMCALAcceptanceTPC( fRecAnalysis->GetCorrectedTotEtEMCALAcceptanceTPC() );
 191         fMCAnalysis->FillSimTotEtMinusRecoTotEtEMCALAcceptanceITS( fRecAnalysis->GetCorrectedTotEtEMCALAcceptanceITS() );
 192         fMCAnalysis->FillSimTotEtMinusRecoTotEtEMCALAcceptanceTPCNoPID( fRecAnalysis->GetCorrectedTotEtEMCALAcceptanceTPCNoPID() );
 193         fMCAnalysis->FillSimTotEtMinusRecoTotEtEMCALAcceptanceITSNoPID( fRecAnalysis->GetCorrectedTotEtEMCALAcceptanceITSNoPID() );
 194         fMCAnalysis->FillSimHadEtMinusRecoHadEtEMCALAcceptanceTPC( fRecAnalysis->GetCorrectedHadEtEMCALAcceptanceTPC() );
 195         fMCAnalysis->FillSimHadEtMinusRecoHadEtEMCALAcceptanceITS( fRecAnalysis->GetCorrectedHadEtEMCALAcceptanceITS() );
 196         fMCAnalysis->FillSimHadEtMinusRecoHadEtEMCALAcceptanceTPCNoPID( fRecAnalysis->GetCorrectedHadEtEMCALAcceptanceTPCNoPID() );
 197         fMCAnalysis->FillSimHadEtMinusRecoHadEtEMCALAcceptanceITSNoPID( fRecAnalysis->GetCorrectedHadEtEMCALAcceptanceITSNoPID() );
 198       }
 199       if(fMCAnalysis->PHOS()){
 200         fMCAnalysis->FillSimTotEtMinusRecoTotEtPHOSAcceptanceTPC( fRecAnalysis->GetCorrectedTotEtPHOSAcceptanceTPC() );
 201         fMCAnalysis->FillSimTotEtMinusRecoTotEtPHOSAcceptanceITS( fRecAnalysis->GetCorrectedTotEtPHOSAcceptanceITS() );
 202         fMCAnalysis->FillSimTotEtMinusRecoTotEtPHOSAcceptanceTPCNoPID( fRecAnalysis->GetCorrectedTotEtPHOSAcceptanceTPCNoPID() );
 203         fMCAnalysis->FillSimTotEtMinusRecoTotEtPHOSAcceptanceITSNoPID( fRecAnalysis->GetCorrectedTotEtPHOSAcceptanceITSNoPID() );
 204         fMCAnalysis->FillSimHadEtMinusRecoHadEtPHOSAcceptanceTPC( fRecAnalysis->GetCorrectedHadEtPHOSAcceptanceTPC() );
 205         fMCAnalysis->FillSimHadEtMinusRecoHadEtPHOSAcceptanceITS( fRecAnalysis->GetCorrectedHadEtPHOSAcceptanceITS() );
 206         fMCAnalysis->FillSimHadEtMinusRecoHadEtPHOSAcceptanceTPCNoPID( fRecAnalysis->GetCorrectedHadEtPHOSAcceptanceTPCNoPID() );
 207         fMCAnalysis->FillSimHadEtMinusRecoHadEtPHOSAcceptanceITSNoPID( fRecAnalysis->GetCorrectedHadEtPHOSAcceptanceITSNoPID() );
 208       }
 209       if(fMCAnalysis->PiKP() && fMCAnalysis->Full()){
 210         fMCAnalysis->FillSimPiKPMinusRecoPiKPFullAcceptanceTPC(fRecAnalysis->GetCorrectedPiKPEtFullAcceptanceTPC());
 211         fMCAnalysis->FillSimPiKPMinusRecoPiKPFullAcceptanceITS(fRecAnalysis->GetCorrectedPiKPEtFullAcceptanceITS());
 212         fMCAnalysis->FillSimPiKPMinusRecoPiKPFullAcceptanceTPCNoPID(fRecAnalysis->GetCorrectedPiKPEtFullAcceptanceTPCNoPID());
 213         fMCAnalysis->FillSimPiKPMinusRecoPiKPFullAcceptanceITSNoPID(fRecAnalysis->GetCorrectedPiKPEtFullAcceptanceITSNoPID());
 214       }
 215     }
 216   }
 217 //cout<<"End Event"<<endl<<endl;
 218 // Post output data.
 219 PostData(1, fOutputList);
 220 }
 221
 222 //________________________________________________________________________
 223 void AliAnalysisTaskHadEt::Terminate(Option_t *)
 224 {
 225     // Draw result to the screen
 226     // Called once at the end of the query
 227
 228     fOutputList = dynamic_cast<TList*> (GetOutputData(1));
 229     if (!fOutputList) {
 230         printf("ERROR: Output list not available\n");
 231         return;
 232     }
 233 }
 234
 235