PWG4/totEt/AliAnalysisHadEtReconstructed.cxx

   1 //_________________________________________________________________________
   2 //  Utility Class for transverse energy studies, charged hadrons
   3 //  Base class for ESD analysis
   4 //  - reconstruction output
   5 // implementation file
   6 //
   7 //Created by Christine Nattrass, Rebecca Scott, Irakli Martashvili
   8 //University of Tennessee at Knoxville
   9 //_________________________________________________________________________
  10
  11 #include <TROOT.h>
  12 #include <TSystem.h>
  13 #include <TInterpreter.h>
  14 #include "AliAnalysisHadEtReconstructed.h"
  15 #include "AliAnalysisEtCuts.h"
  16 #include "AliESDtrack.h"
  17 #include "AliESDCaloCluster.h"
  18 #include "AliVEvent.h"
  19 #include "AliESDEvent.h"
  20 #include "AliESDtrackCuts.h"
  21 #include "AliESDpid.h"
  22 #include "AliVParticle.h"
  23 #include <iostream>
  24 #include "AliAnalysisHadEtCorrections.h"
  25 #include "TFile.h"
  26 #include "TString.h"
  27
  28 using namespace std;
  29
  30 ClassImp(AliAnalysisHadEtReconstructed);
  31
  32
  33 AliAnalysisHadEtReconstructed::AliAnalysisHadEtReconstructed() :
  34         AliAnalysisHadEt()
  35         ,fCorrections(0)
  36         ,fConfigFile("ConfigHadEtAnalysis.C")
  37     ,fCorrTotEtFullAcceptanceTPC(0)
  38     ,fCorrTotEtFullAcceptanceITS(0)
  39     ,fCorrHadEtFullAcceptanceTPC(0)
  40     ,fCorrHadEtFullAcceptanceITS(0)
  41     ,fCorrTotEtEMCALAcceptanceTPC(0)
  42     ,fCorrTotEtEMCALAcceptanceITS(0)
  43     ,fCorrHadEtEMCALAcceptanceTPC(0)
  44     ,fCorrHadEtEMCALAcceptanceITS(0)
  45     ,fCorrTotEtPHOSAcceptanceTPC(0)
  46     ,fCorrTotEtPHOSAcceptanceITS(0)
  47     ,fCorrHadEtPHOSAcceptanceTPC(0)
  48     ,fCorrHadEtPHOSAcceptanceITS(0)
  49     ,fCorrectedHadEtFullAcceptanceTPCNoPID(0)
  50     ,fCorrectedHadEtFullAcceptanceITSNoPID(0)
  51     ,fCorrectedHadEtEMCALAcceptanceTPCNoPID(0)
  52     ,fCorrectedHadEtEMCALAcceptanceITSNoPID(0)
  53     ,fCorrectedHadEtPHOSAcceptanceTPCNoPID(0)
  54     ,fCorrectedHadEtPHOSAcceptanceITSNoPID(0)
  55     ,fCorrectedHadEtFullAcceptanceTPC(0)
  56     ,fCorrectedHadEtFullAcceptanceITS(0)
  57     ,fCorrectedHadEtEMCALAcceptanceTPC(0)
  58     ,fCorrectedHadEtEMCALAcceptanceITS(0)
  59     ,fCorrectedHadEtPHOSAcceptanceTPC(0)
  60     ,fCorrectedHadEtPHOSAcceptanceITS(0)
  61     ,fRawEtFullAcceptanceTPC(0)
  62     ,fRawEtFullAcceptanceITS(0)
  63     ,fRawEtEMCALAcceptanceTPC(0)
  64     ,fRawEtEMCALAcceptanceITS(0)
  65     ,fRawEtPHOSAcceptanceTPC(0)
  66     ,fRawEtPHOSAcceptanceITS(0)
  67     ,fRawEtFullAcceptanceTPCNoPID(0)
  68     ,fRawEtFullAcceptanceITSNoPID(0)
  69     ,fRawEtEMCALAcceptanceTPCNoPID(0)
  70     ,fRawEtEMCALAcceptanceITSNoPID(0)
  71     ,fRawEtPHOSAcceptanceTPCNoPID(0)
  72     ,fRawEtPHOSAcceptanceITSNoPID(0)
  73 {
  74
  75 }
  76
  77 AliAnalysisHadEtReconstructed::~AliAnalysisHadEtReconstructed()
  78 {
  79   delete fCorrections;
  80 }
  81
  82 Int_t AliAnalysisHadEtReconstructed::AnalyseEvent(AliVEvent* ev)
  83 { // analyse ESD event
  84     ResetEventValues();
  85
  86     AliESDEvent *realEvent = dynamic_cast<AliESDEvent*>(ev);
  87     //for PID
  88     AliESDpid *pID = new AliESDpid();
  89     pID->MakePID(realEvent);
  90
  91     TString *strTPC = new TString("TPC");
  92     TString *strITS = new TString("ITS");
  93     TString *strTPCITS = new TString("TPCITS");
  94     bool isTPC = false;
  95     for(Int_t cutset=0;cutset<2;cutset++){
  96       TString *cutName;
  97       TObjArray* list;
  98       switch(cutset){
  99       case 0:
 100         cutName = strTPC;
 101         list = fEsdtrackCutsTPC->GetAcceptedTracks(realEvent);
 102         isTPC = true;
 103         break;
 104       case 1:
 105         cutName = strITS;
 106         list = fEsdtrackCutsITS->GetAcceptedTracks(realEvent);
 107         break;
 108       case 2:
 109         cutName = strTPCITS;
 110         list = fEsdtrackCutsITSTPC->GetAcceptedTracks(realEvent);
 111         break;
 112       default:
 113         cerr<<"Error:  cannot fill histograms!"<<endl;
 114         return -1;
 115       }
 116       Int_t nGoodTracks = list->GetEntries();
 117       //cout<<nGoodTracks<<" "<<cutName->Data()<<" tracks"<<endl;
 118       for (Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++)
 119         {
 120
 121
 122           AliESDtrack *track = dynamic_cast<AliESDtrack*> (list->At(iTrack));
 123           if (!track)
 124             {
 125               Printf("ERROR: Could not get track %d", iTrack);
 126               continue;
 127             }
 128           else{
 129             if(TMath::Abs(track->Eta())>fCorrections->GetEtaCut()) continue;
 130             Float_t nSigmaPion,nSigmaProton,nSigmaKaon,nSigmaElectron;
 131
 132             if(cutset!=1){
 133               nSigmaPion = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kPion));
 134               nSigmaProton = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kProton));
 135               nSigmaKaon = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kKaon));
 136               nSigmaElectron = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kElectron));
 137             }
 138             else{
 139               nSigmaPion = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kPion));
 140               nSigmaProton = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kProton));
 141               nSigmaKaon = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kKaon));
 142               nSigmaElectron = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kElectron));
 143             }
 144             bool isPion = (nSigmaPion<3.0 && nSigmaProton>2.0 && nSigmaKaon>2.0);
 145             bool isElectron = (nSigmaElectron<2.0 && nSigmaPion>4.0 && nSigmaProton>3.0 && nSigmaKaon>3.0);
 146             bool isKaon = (nSigmaPion>3.0 && nSigmaProton>2.0 && nSigmaKaon<2.0);
 147             bool isProton = (nSigmaPion>3.0 && nSigmaProton<2.0 && nSigmaKaon>2.0);
 148
 149             //bool IsElectron = false;
 150             bool unidentified = (!isProton && !isKaon && !isElectron);
 151             Float_t dEdx = track->GetTPCsignal();
 152             if(cutset==1) dEdx = track->GetITSsignal();
 153             FillHisto2D(Form("dEdxDataAll%s",cutName->Data()),track->P(),dEdx,1.0);
 154
 155             bool inPHOS = IsInPHOS(track);
 156             bool inEMCAL = IsInEMCAL(track);
 157             //if(!(fCorrections->GetEfficiencyPionTPC())) cerr<<"Uh-oh!  No histogram!"<<endl;
 158
 159             Float_t corrBkgd=0.0;
 160             Float_t corrNotID=0.0;
 161             Float_t corrNoID = fCorrections->GetNotIDCorrectionNoPID(track->Pt());
 162             Float_t corrEff = 0.0;
 163             Float_t corrEffNoID = 0.0;
 164             if(cutset==0){//TPC
 165               corrBkgd = fCorrections->GetBackgroundCorrectionTPC(track->Pt());
 166               corrEffNoID = fCorrections->GetTPCEfficiencyCorrectionHadron(track->Pt());
 167               corrNotID = fCorrections->GetNotIDCorrectionTPC(track->Pt());
 168             }
 169             if(cutset==1){//ITS
 170               corrBkgd = fCorrections->GetBackgroundCorrectionITS(track->Pt());
 171               //corrEffNoID = fCorrections->GetITSEfficiencyCorrectionHadron(track->Pt());
 172               corrNotID = fCorrections->GetNotIDCorrectionITS(track->Pt());
 173             }
 174             Float_t et = 0.0;
 175             Float_t etNoID = Et(track->P(),track->Theta(),fPiPlusCode,track->Charge());
 176             Float_t etpartialcorrected = 0.0;
 177             Float_t etpartialcorrectedNoID = corrNoID*corrBkgd*corrEffNoID*etNoID;
 178             FillHisto2D(Form("EtDataRaw%sNoID",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrectedNoID);
 179
 180             if(isPion){
 181               FillHisto2D(Form("dEdxDataPion%s",cutName->Data()),track->P(),dEdx,1.0);
 182               et = Et(track->P(),track->Theta(),fPiPlusCode,track->Charge());
 183               if(cutset==0){corrEff = fCorrections->GetTPCEfficiencyCorrectionPion(track->Pt());}
 184               //else{corrEff = fCorrections->GetITSEfficiencyCorrectionPion(track->Pt());}
 185               etpartialcorrected = et*corrBkgd*corrEff;
 186
 187               if(track->Charge()>0.0){
 188                 FillHisto2D(Form("EtDataRaw%sPiPlus",cutName->Data()),track->Pt(),track->Eta(),et);
 189                 FillHisto2D(Form("EtDataCorrected%sPiPlus",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
 190               }
 191               else{
 192                 FillHisto2D(Form("EtDataRaw%sPiMinus",cutName->Data()),track->Pt(),track->Eta(),et);
 193                 FillHisto2D(Form("EtDataCorrected%sPiMinus",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
 194               }
 195             }
 196             if(isKaon){
 197               FillHisto2D(Form("dEdxDataKaon%s",cutName->Data()),track->P(),dEdx,1.0);
 198               et = Et(track->P(),track->Theta(),fKPlusCode,track->Charge());
 199               if(cutset==0){corrEff = fCorrections->GetTPCEfficiencyCorrectionKaon(track->Pt());}
 200               //else{corrEff = fCorrections->GetITSEfficiencyCorrectionKaon(track->Pt());}
 201               etpartialcorrected = et*corrBkgd*corrEff;
 202
 203               if(track->Charge()>0.0){
 204                 FillHisto2D(Form("EtDataRaw%sKPlus",cutName->Data()),track->Pt(),track->Eta(),et);
 205                 FillHisto2D(Form("EtDataCorrected%sKPlus",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
 206               }
 207               else{
 208                 FillHisto2D(Form("EtDataRaw%sKMinus",cutName->Data()),track->Pt(),track->Eta(),et);
 209                 FillHisto2D(Form("EtDataCorrected%sKMinus",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
 210               }
 211             }
 212             if(isProton){
 213               FillHisto2D(Form("dEdxDataProton%s",cutName->Data()),track->P(),dEdx,1.0);
 214               et = Et(track->P(),track->Theta(),fProtonCode,track->Charge());
 215               if(cutset==0){corrEff = fCorrections->GetTPCEfficiencyCorrectionProton(track->Pt());}
 216               //else{corrEff = fCorrections->GetITSEfficiencyCorrectionProton(track->Pt());}
 217               etpartialcorrected = et*corrBkgd*corrEff;
 218
 219               if(track->Charge()>0.0){
 220                 FillHisto2D(Form("EtDataRaw%sProton",cutName->Data()),track->Pt(),track->Eta(),et);
 221                 FillHisto2D(Form("EtDataCorrected%sProton",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
 222               }
 223               else{
 224                 FillHisto2D(Form("EtDataRaw%sAntiProton",cutName->Data()),track->Pt(),track->Eta(),et);
 225                 FillHisto2D(Form("EtDataCorrected%sAntiProton",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
 226               }
 227             }
 228             if(isElectron){
 229               FillHisto2D(Form("dEdxDataProton%s",cutName->Data()),track->P(),dEdx,1.0);
 230               //et = Et(track->P(),track->Theta(),fPiPlusCode,track->Charge());
 231             }
 232             if(unidentified){
 233               FillHisto2D(Form("dEdxDataUnidentified%s",cutName->Data()),track->P(),dEdx,1.0);
 234               et = Et(track->P(),track->Theta(),fPiPlusCode,track->Charge());
 235               etpartialcorrected = et*corrBkgd*corrEffNoID*corrNotID;
 236               FillHisto2D(Form("EtDataCorrected%sUnidentified",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
 237             }
 238             if(!isTPC) etpartialcorrected = etpartialcorrectedNoID;//Not using PID for ITS
 239             AddEt(et,etNoID,etpartialcorrected,etpartialcorrectedNoID,track->Pt(),isTPC,inPHOS,inEMCAL);
 240             //if(inEMCAL) cout<<"I should add a track"<<endl;
 241           }
 242         }
 243     }
 244 //   cout<<"Finishing with Raw/Corrected Et in full, PHOS, EMCAL acceptance of "
 245 //       << GetRawEtFullAcceptanceITS() <<"/"
 246 //       << GetCorrectedHadEtFullAcceptanceITS() <<", "
 247 //       << GetRawEtPHOSAcceptanceITS() <<"/"
 248 //       << GetCorrectedHadEtPHOSAcceptanceITS() <<", "
 249 //       << GetRawEtEMCALAcceptanceITS() <<"/"
 250 //       << GetCorrectedHadEtEMCALAcceptanceITS() <<endl;
 251 //   cout<<"Finishing with Raw/Corrected Et w/o PID in full, PHOS, EMCAL acceptance of "
 252 //       << GetRawEtFullAcceptanceITSNoPID() <<"/"
 253 //       << GetCorrectedHadEtFullAcceptanceITSNoPID() <<", "
 254 //       << GetRawEtPHOSAcceptanceITSNoPID() <<"/"
 255 //       << GetCorrectedHadEtPHOSAcceptanceITSNoPID() <<", "
 256 //       << GetRawEtEMCALAcceptanceITSNoPID() <<"/"
 257 //       << GetCorrectedHadEtEMCALAcceptanceITSNoPID() <<endl;
 258 //   cout<<"Finishing with Raw/Corrected Et in full, PHOS, EMCAL acceptance of "
 259 //       << GetRawEtFullAcceptanceTPC() <<"/"
 260 //       << GetCorrectedHadEtFullAcceptanceTPC() <<", "
 261 //       << GetRawEtPHOSAcceptanceTPC() <<"/"
 262 //       << GetCorrectedHadEtPHOSAcceptanceTPC() <<", "
 263 //       << GetRawEtEMCALAcceptanceTPC() <<"/"
 264 //       << GetCorrectedHadEtEMCALAcceptanceTPC() <<endl;
 265 //   cout<<"Finishing with Raw/Corrected Et w/o PID in full, PHOS, EMCAL acceptance of "
 266 //       << GetRawEtFullAcceptanceTPCNoPID() <<"/"
 267 //       << GetCorrectedHadEtFullAcceptanceTPCNoPID() <<", "
 268 //       << GetRawEtPHOSAcceptanceTPCNoPID() <<"/"
 269 //       << GetCorrectedHadEtPHOSAcceptanceTPCNoPID() <<", "
 270 //       << GetRawEtEMCALAcceptanceTPCNoPID() <<"/"
 271 //       << GetCorrectedHadEtEMCALAcceptanceTPCNoPID() <<endl;
 272 //   cout<<"Correction factors "
 273 //       <<fCorrTotEtFullAcceptanceTPC<<", "<<fCorrTotEtFullAcceptanceITS<<", "<<fCorrHadEtFullAcceptanceTPC<<", "<<fCorrHadEtFullAcceptanceITS<<","
 274 //       <<fCorrTotEtEMCALAcceptanceTPC<<", "<<fCorrTotEtEMCALAcceptanceITS<<", "<<fCorrHadEtEMCALAcceptanceTPC<<", "<<fCorrHadEtEMCALAcceptanceITS<<","
 275 //       <<fCorrTotEtPHOSAcceptanceTPC<<", "<<fCorrTotEtPHOSAcceptanceITS<<", "<<fCorrHadEtPHOSAcceptanceTPC<<", "<<fCorrHadEtPHOSAcceptanceITS<<endl;
 276     FillHisto1D("RecoHadEtFullAcceptanceTPC",GetCorrectedHadEtFullAcceptanceTPC(),1.0);
 277     FillHisto1D("RecoTotEtFullAcceptanceTPC",GetCorrectedTotEtFullAcceptanceTPC(),1.0);
 278     FillHisto1D("RecoHadEtEMCALAcceptanceTPC",GetCorrectedHadEtEMCALAcceptanceTPC(),1.0);
 279     FillHisto1D("RecoTotEtEMCALAcceptanceTPC",GetCorrectedTotEtEMCALAcceptanceTPC(),1.0);
 280     FillHisto1D("RecoHadEtPHOSAcceptanceTPC",GetCorrectedHadEtPHOSAcceptanceTPC(),1.0);
 281     FillHisto1D("RecoTotEtPHOSAcceptanceTPC",GetCorrectedTotEtPHOSAcceptanceTPC(),1.0);
 282     FillHisto1D("RecoHadEtFullAcceptanceTPCNoPID",GetCorrectedHadEtFullAcceptanceTPCNoPID(),1.0);
 283     FillHisto1D("RecoTotEtFullAcceptanceTPCNoPID",GetCorrectedTotEtFullAcceptanceTPCNoPID(),1.0);
 284     FillHisto1D("RecoHadEtEMCALAcceptanceTPCNoPID",GetCorrectedHadEtEMCALAcceptanceTPCNoPID(),1.0);
 285     FillHisto1D("RecoTotEtEMCALAcceptanceTPCNoPID",GetCorrectedTotEtEMCALAcceptanceTPCNoPID(),1.0);
 286     FillHisto1D("RecoHadEtPHOSAcceptanceTPCNoPID",GetCorrectedHadEtPHOSAcceptanceTPCNoPID(),1.0);
 287     FillHisto1D("RecoTotEtPHOSAcceptanceTPCNoPID",GetCorrectedTotEtPHOSAcceptanceTPCNoPID(),1.0);
 288     FillHisto1D("RecoHadEtFullAcceptanceITS",GetCorrectedHadEtFullAcceptanceITS(),1.0);
 289     FillHisto1D("RecoTotEtFullAcceptanceITS",GetCorrectedTotEtFullAcceptanceITS(),1.0);
 290     FillHisto1D("RecoHadEtEMCALAcceptanceITS",GetCorrectedHadEtEMCALAcceptanceITS(),1.0);
 291     FillHisto1D("RecoTotEtEMCALAcceptanceITS",GetCorrectedTotEtEMCALAcceptanceITS(),1.0);
 292     FillHisto1D("RecoHadEtPHOSAcceptanceITS",GetCorrectedHadEtPHOSAcceptanceITS(),1.0);
 293     FillHisto1D("RecoTotEtPHOSAcceptanceITS",GetCorrectedTotEtPHOSAcceptanceITS(),1.0);
 294     FillHisto1D("RecoHadEtFullAcceptanceITSNoPID",GetCorrectedHadEtFullAcceptanceITSNoPID(),1.0);
 295     FillHisto1D("RecoTotEtFullAcceptanceITSNoPID",GetCorrectedTotEtFullAcceptanceITSNoPID(),1.0);
 296     FillHisto1D("RecoHadEtEMCALAcceptanceITSNoPID",GetCorrectedHadEtEMCALAcceptanceITSNoPID(),1.0);
 297     FillHisto1D("RecoTotEtEMCALAcceptanceITSNoPID",GetCorrectedTotEtEMCALAcceptanceITSNoPID(),1.0);
 298     FillHisto1D("RecoHadEtPHOSAcceptanceITSNoPID",GetCorrectedHadEtPHOSAcceptanceITSNoPID(),1.0);
 299     FillHisto1D("RecoTotEtPHOSAcceptanceITSNoPID",GetCorrectedTotEtPHOSAcceptanceITSNoPID(),1.0);
 300
 301     FillHisto1D("RecoRawEtFullAcceptanceTPC",GetRawEtFullAcceptanceTPC(),1.0);
 302     FillHisto1D("RecoRawEtEMCALAcceptanceTPC",GetRawEtEMCALAcceptanceTPC(),1.0);
 303     FillHisto1D("RecoRawEtPHOSAcceptanceTPC",GetRawEtPHOSAcceptanceTPC(),1.0);
 304     FillHisto1D("RecoRawEtFullAcceptanceTPCNoPID",GetRawEtFullAcceptanceTPCNoPID(),1.0);
 305     FillHisto1D("RecoRawEtEMCALAcceptanceTPCNoPID",GetRawEtEMCALAcceptanceTPCNoPID(),1.0);
 306     FillHisto1D("RecoRawEtPHOSAcceptanceTPCNoPID",GetRawEtPHOSAcceptanceTPCNoPID(),1.0);
 307     FillHisto1D("RecoRawEtFullAcceptanceITS",GetRawEtFullAcceptanceITS(),1.0);
 308     FillHisto1D("RecoRawEtEMCALAcceptanceITS",GetRawEtEMCALAcceptanceITS(),1.0);
 309     FillHisto1D("RecoRawEtPHOSAcceptanceITS",GetRawEtPHOSAcceptanceITS(),1.0);
 310     FillHisto1D("RecoRawEtFullAcceptanceITSNoPID",GetRawEtFullAcceptanceITSNoPID(),1.0);
 311     FillHisto1D("RecoRawEtEMCALAcceptanceITSNoPID",GetRawEtEMCALAcceptanceITSNoPID(),1.0);
 312     FillHisto1D("RecoRawEtPHOSAcceptanceITSNoPID",GetRawEtPHOSAcceptanceITSNoPID(),1.0);
 313     delete pID;
 314     delete strTPC;
 315     delete strITS;
 316     delete strTPCITS;
 317     return 1;
 318 }
 319 void AliAnalysisHadEtReconstructed::AddEt(Float_t rawEt, Float_t rawEtNoPID, Float_t corrEt, Float_t corrEtNoPID, Float_t pt, Bool_t IsTPC, Bool_t InPHOS, Bool_t InEMCAL) {//Adding Et to each of the variables that tracks et event by event
 320   if(pt>=AliAnalysisHadEt::fgPtTPCCutOff && IsTPC){//TPC tracks
 321     //adding to the raw Et
 322     //if(InEMCAL) cout<<"Adding "<<rawEt<<" to the raw Et"<<endl;
 323     fRawEtFullAcceptanceTPC += rawEt;
 324     if(InPHOS)fRawEtPHOSAcceptanceTPC += rawEt;
 325     if(InEMCAL)fRawEtEMCALAcceptanceTPC += rawEt;
 326     fRawEtFullAcceptanceTPCNoPID += rawEtNoPID;
 327     if(InPHOS)fRawEtPHOSAcceptanceTPCNoPID += rawEtNoPID;
 328     if(InEMCAL)fRawEtEMCALAcceptanceTPCNoPID += rawEtNoPID;
 329     //adding to the corrected Et
 330     //if(InPHOS) cout<<"Adding "<<corrEt<<" to the corrected Et"<<endl;
 331     fCorrectedHadEtFullAcceptanceTPC += corrEt;
 332     if(InPHOS)fCorrectedHadEtPHOSAcceptanceTPC += corrEt;
 333     if(InEMCAL)fCorrectedHadEtEMCALAcceptanceTPC += corrEt;
 334     fCorrectedHadEtFullAcceptanceTPCNoPID += corrEtNoPID;
 335     if(InPHOS)fCorrectedHadEtPHOSAcceptanceTPCNoPID += corrEtNoPID;
 336     if(InEMCAL)fCorrectedHadEtEMCALAcceptanceTPCNoPID += corrEtNoPID;
 337   }
 338   if(pt<AliAnalysisHadEt::fgPtTPCCutOff &&pt>=AliAnalysisHadEt::fgPtITSCutOff && !IsTPC){//ITS tracks
 339     //adding to the raw Et
 340     fRawEtFullAcceptanceITS += rawEt;
 341     if(InPHOS)fRawEtPHOSAcceptanceITS += rawEt;
 342     if(InEMCAL)fRawEtEMCALAcceptanceITS += rawEt;
 343     fRawEtFullAcceptanceITSNoPID += rawEtNoPID;
 344     if(InPHOS)fRawEtPHOSAcceptanceITSNoPID += rawEtNoPID;
 345     if(InEMCAL)fRawEtEMCALAcceptanceITSNoPID += rawEtNoPID;
 346     //adding to the corrected Et
 347     fCorrectedHadEtFullAcceptanceITS += corrEt;
 348     if(InPHOS)fCorrectedHadEtPHOSAcceptanceITS += corrEt;
 349     if(InEMCAL)fCorrectedHadEtEMCALAcceptanceITS += corrEt;
 350     fCorrectedHadEtFullAcceptanceITSNoPID += corrEtNoPID;
 351     if(InPHOS)fCorrectedHadEtPHOSAcceptanceITSNoPID += corrEtNoPID;
 352     if(InEMCAL)fCorrectedHadEtEMCALAcceptanceITSNoPID += corrEtNoPID;
 353   }
 354 }
 355
 356 Bool_t AliAnalysisHadEtReconstructed::IsInPHOS(AliESDtrack *track){//This function will need to be elaborated on later to include PHOS dead channels
 357   return   TMath::Abs(track->Eta()) < fCuts->GetGeometryPhosEtaAccCut()//in eta acceptance
 358     && track->Phi()*180.0/TMath::Pi() > fCuts->GetGeometryPhosPhiAccMinCut()//greater than the minimum phi
 359     && track->Phi()*180.0/TMath::Pi() < fCuts->GetGeometryPhosPhiAccMaxCut();//less than the maximum phi
 360 }
 361 Bool_t AliAnalysisHadEtReconstructed::IsInEMCAL(AliESDtrack *track){//This function will need to be elaborated on later to include EMCAL dead channels
 362   //cout<<"Eta: |"<<track->Eta()<<"|<"<< fCuts->GetGeometryEmcalEtaAccCut() <<"; phi: "<<fCuts->GetGeometryEmcalPhiAccMinCut()<<"<"<<track->Phi()*180.0/TMath::Pi()<<"<"<<fCuts->GetGeometryEmcalPhiAccMaxCut()<<endl;
 363   return   TMath::Abs(track->Eta()) < fCuts->GetGeometryEmcalEtaAccCut()//in eta acceptance
 364     && track->Phi()*180.0/TMath::Pi() > fCuts->GetGeometryEmcalPhiAccMinCut()//greater than the minimum phi
 365     && track->Phi()*180.0/TMath::Pi() < fCuts->GetGeometryEmcalPhiAccMaxCut();//less than the maximum phi
 366 }
 367 Bool_t AliAnalysisHadEtReconstructed::CheckGoodVertex(AliVParticle* track)
 368 { // check vertex
 369
 370     Float_t bxy = 999.;
 371     Float_t bz = 999.;
 372     dynamic_cast<AliESDtrack*>(track)->GetImpactParametersTPC(bxy,bz);
 373
 374     bool status = (TMath::Abs(track->Xv()) < fCuts->GetReconstructedVertexXCut()) &&
 375       (TMath::Abs(track->Yv()) < fCuts->GetReconstructedVertexYCut()) &&
 376       (TMath::Abs(track->Zv()) < fCuts->GetReconstructedVertexZCut()) &&
 377       (TMath::Abs(bxy) < fCuts->GetReconstructedIPxyCut()) &&
 378       (TMath::Abs(bz) < fCuts->GetReconstructedIPzCut());
 379
 380     return status;
 381 }
 382
 383 void AliAnalysisHadEtReconstructed::Init()
 384 { // Init
 385     AliAnalysisHadEt::Init();
 386
 387   if (fConfigFile.Length()) {
 388     gROOT->LoadMacro(fConfigFile);
 389     fCorrections = (AliAnalysisHadEtCorrections *) gInterpreter->ProcessLine("ConfigHadEtAnalysis()");
 390     fCorrTotEtFullAcceptanceTPC = fCorrections->GetConstantCorrections(kTRUE,fgPtTPCCutOff,"Full");
 391     fCorrTotEtFullAcceptanceITS = fCorrections->GetConstantCorrections(kTRUE,fgPtITSCutOff,"Full");
 392     fCorrHadEtFullAcceptanceTPC = fCorrections->GetConstantCorrections(kFALSE,fgPtTPCCutOff,"Full");
 393     fCorrHadEtFullAcceptanceITS = fCorrections->GetConstantCorrections(kFALSE,fgPtITSCutOff,"Full");
 394     fCorrTotEtEMCALAcceptanceTPC = fCorrections->GetConstantCorrections(kTRUE,fgPtTPCCutOff,"EMCAL");
 395     fCorrTotEtEMCALAcceptanceITS = fCorrections->GetConstantCorrections(kTRUE,fgPtITSCutOff,"EMCAL");
 396     fCorrHadEtEMCALAcceptanceTPC = fCorrections->GetConstantCorrections(kFALSE,fgPtTPCCutOff,"EMCAL");
 397     fCorrHadEtEMCALAcceptanceITS = fCorrections->GetConstantCorrections(kFALSE,fgPtITSCutOff,"EMCAL");
 398     fCorrTotEtPHOSAcceptanceTPC = fCorrections->GetConstantCorrections(kTRUE,fgPtTPCCutOff,"PHOS");
 399     fCorrTotEtPHOSAcceptanceITS = fCorrections->GetConstantCorrections(kTRUE,fgPtITSCutOff,"PHOS");
 400     fCorrHadEtPHOSAcceptanceTPC = fCorrections->GetConstantCorrections(kFALSE,fgPtTPCCutOff,"PHOS");
 401     fCorrHadEtPHOSAcceptanceITS = fCorrections->GetConstantCorrections(kFALSE,fgPtITSCutOff,"PHOS");
 402
 403   }
 404 }
 405
 406 void AliAnalysisHadEtReconstructed::ResetEventValues(){//resetting event by event et's
 407   AliAnalysisHadEt::ResetEventValues();
 408      fCorrectedHadEtFullAcceptanceTPCNoPID=0.0;
 409      fCorrectedHadEtFullAcceptanceITSNoPID=0.0;
 410      fCorrectedHadEtEMCALAcceptanceTPCNoPID=0.0;
 411      fCorrectedHadEtEMCALAcceptanceITSNoPID=0.0;
 412      fCorrectedHadEtPHOSAcceptanceTPCNoPID=0.0;
 413      fCorrectedHadEtPHOSAcceptanceITSNoPID=0.0;
 414      fCorrectedHadEtFullAcceptanceTPC=0.0;
 415      fCorrectedHadEtFullAcceptanceITS=0.0;
 416      fCorrectedHadEtEMCALAcceptanceTPC=0.0;
 417      fCorrectedHadEtEMCALAcceptanceITS=0.0;
 418      fCorrectedHadEtPHOSAcceptanceTPC=0.0;
 419      fCorrectedHadEtPHOSAcceptanceITS=0.0;
 420      fRawEtFullAcceptanceTPC=0.0;
 421      fRawEtFullAcceptanceITS=0.0;
 422      fRawEtEMCALAcceptanceTPC=0.0;
 423      fRawEtEMCALAcceptanceITS=0.0;
 424      fRawEtPHOSAcceptanceTPC=0.0;
 425      fRawEtPHOSAcceptanceITS=0.0;
 426      fRawEtFullAcceptanceTPCNoPID=0.0;
 427      fRawEtFullAcceptanceITSNoPID=0.0;
 428      fRawEtEMCALAcceptanceTPCNoPID=0.0;
 429      fRawEtEMCALAcceptanceITSNoPID=0.0;
 430      fRawEtPHOSAcceptanceTPCNoPID=0.0;
 431      fRawEtPHOSAcceptanceITSNoPID=0.0;
 432
 433      if(TMath::Abs(fCorrTotEtFullAcceptanceTPC)<1e-3){
 434        if (fConfigFile.Length()) {
 435          cout<<"Rereading fCorrections file..."<<endl;
 436          gROOT->LoadMacro(fConfigFile);
 437          fCorrections = (AliAnalysisHadEtCorrections *) gInterpreter->ProcessLine("ConfigHadEtAnalysis()");
 438          fCorrTotEtFullAcceptanceTPC = fCorrections->GetConstantCorrections(kTRUE,fgPtTPCCutOff,"Full");
 439          fCorrTotEtFullAcceptanceITS = fCorrections->GetConstantCorrections(kTRUE,fgPtITSCutOff,"Full");
 440          fCorrHadEtFullAcceptanceTPC = fCorrections->GetConstantCorrections(kFALSE,fgPtTPCCutOff,"Full");
 441          fCorrHadEtFullAcceptanceITS = fCorrections->GetConstantCorrections(kFALSE,fgPtITSCutOff,"Full");
 442          fCorrTotEtEMCALAcceptanceTPC = fCorrections->GetConstantCorrections(kTRUE,fgPtTPCCutOff,"EMCAL");
 443          fCorrTotEtEMCALAcceptanceITS = fCorrections->GetConstantCorrections(kTRUE,fgPtITSCutOff,"EMCAL");
 444          fCorrHadEtEMCALAcceptanceTPC = fCorrections->GetConstantCorrections(kFALSE,fgPtTPCCutOff,"EMCAL");
 445          fCorrHadEtEMCALAcceptanceITS = fCorrections->GetConstantCorrections(kFALSE,fgPtITSCutOff,"EMCAL");
 446          fCorrTotEtPHOSAcceptanceTPC = fCorrections->GetConstantCorrections(kTRUE,fgPtTPCCutOff,"PHOS");
 447          fCorrTotEtPHOSAcceptanceITS = fCorrections->GetConstantCorrections(kTRUE,fgPtITSCutOff,"PHOS");
 448          fCorrHadEtPHOSAcceptanceTPC = fCorrections->GetConstantCorrections(kFALSE,fgPtTPCCutOff,"PHOS");
 449          fCorrHadEtPHOSAcceptanceITS = fCorrections->GetConstantCorrections(kFALSE,fgPtITSCutOff,"PHOS");
 450        }
 451        else{cerr<<"Uh-oh!  Unable to open configuration file!"<<endl;}
 452      }
 453
 454 }
 455 void AliAnalysisHadEtReconstructed::CreateHistograms(){//Creating histograms and adding them to the output TList
 456
 457   TString *strTPC = new TString("TPC");
 458   TString *strITS = new TString("ITS");
 459   TString *strTPCITS = new TString("TPCITS");
 460   for(Int_t i=0;i<2;i++){
 461     TString *cutName;
 462     Float_t maxPtdEdx = 10;
 463     Float_t mindEdx = 35;
 464     Float_t maxdEdx = 150.0;
 465     switch(i){
 466     case 0:
 467       cutName = strTPC;
 468       break;
 469     case 1:
 470       cutName = strITS;
 471       maxPtdEdx = 5;
 472       maxdEdx = 500.0;
 473       break;
 474     case 2:
 475       cutName = strTPCITS;
 476       break;
 477     default:
 478       cerr<<"Error:  cannot make histograms!"<<endl;
 479       return;
 480     }
 481
 482     CreateEtaPtHisto2D(Form("EtDataRaw%sPiPlus",cutName->Data()),"Raw reconstructed E_{T} from identified #pi^{+}");
 483     CreateEtaPtHisto2D(Form("EtDataRaw%sPiMinus",cutName->Data()),"Raw reconstructed E_{T} from identified #pi^{-}");
 484     CreateEtaPtHisto2D(Form("EtDataRaw%sKPlus",cutName->Data()),"Raw reconstructed E_{T} from identified K^{+}");
 485 //     CreateEtaPtHisto2D(Form("EtDataRaw%sEMinus",cutName->Data()),"Raw reconstructed E_{T} from identified e^{-}");
 486 //     CreateEtaPtHisto2D(Form("EtDataRaw%sEPlus",cutName->Data()),"Raw reconstructed E_{T} from identified e^{+}");
 487     CreateEtaPtHisto2D(Form("EtDataRaw%sKMinus",cutName->Data()),"Raw reconstructed E_{T} from identified K^{-}");
 488     CreateEtaPtHisto2D(Form("EtDataRaw%sProton",cutName->Data()),"Raw reconstructed E_{T} from identified p");
 489     CreateEtaPtHisto2D(Form("EtDataRaw%sAntiProton",cutName->Data()),"Raw reconstructed E_{T} from identified #bar{p}");
 490     CreateEtaPtHisto2D(Form("EtDataRaw%sUnidentified",cutName->Data()),"Raw reconstructed E_{T} from unidentified particles using real mass");
 491     CreateEtaPtHisto2D(Form("EtDataRaw%sNoID",cutName->Data()),"Raw reconstructed E_{T} from unidentified particles using real mass");
 492
 493     CreateEtaPtHisto2D(Form("EtDataCorrected%sPiPlus",cutName->Data()),"Corrected reconstructed E_{T} from identified #pi^{+}");
 494     CreateEtaPtHisto2D(Form("EtDataCorrected%sPiMinus",cutName->Data()),"Corrected reconstructed E_{T} from identified #pi^{-}");
 495     CreateEtaPtHisto2D(Form("EtDataCorrected%sKPlus",cutName->Data()),"Corrected reconstructed E_{T} from identified K^{+}");
 496 //     CreateEtaPtHisto2D(Form("EtDataCorrected%sEMinus",cutName->Data()),"Corrected reconstructed E_{T} from identified e^{-}");
 497 //     CreateEtaPtHisto2D(Form("EtDataCorrected%sEPlus",cutName->Data()),"Corrected reconstructed E_{T} from identified e^{+}");
 498     CreateEtaPtHisto2D(Form("EtDataCorrected%sKMinus",cutName->Data()),"Corrected reconstructed E_{T} from identified K^{-}");
 499     CreateEtaPtHisto2D(Form("EtDataCorrected%sProton",cutName->Data()),"Corrected reconstructed E_{T} from identified p");
 500     CreateEtaPtHisto2D(Form("EtDataCorrected%sAntiProton",cutName->Data()),"Corrected reconstructed E_{T} from identified #bar{p}");
 501     CreateEtaPtHisto2D(Form("EtDataCorrected%sUnidentified",cutName->Data()),"Corrected reconstructed E_{T} from unidentified particles using real mass");
 502     CreateEtaPtHisto2D(Form("EtDataCorrected%sNoID",cutName->Data()),"Corrected reconstructed E_{T} from unidentified particles using real mass");
 503
 504
 505     CreateEtaPtHisto2D(Form("EtNData%sPiPlus",cutName->Data()),"Number of reconstructed #pi^{+}");
 506     CreateEtaPtHisto2D(Form("EtNData%sPiMinus",cutName->Data()),"Number of reconstructed #pi^{-}");
 507     CreateEtaPtHisto2D(Form("EtNData%sKPlus",cutName->Data()),"Number of reconstructed K^{+}");
 508     CreateEtaPtHisto2D(Form("EtNData%sKMinus",cutName->Data()),"Number of reconstructed K^{-}");
 509     CreateEtaPtHisto2D(Form("EtNData%sProton",cutName->Data()),"Number of reconstructed p");
 510     CreateEtaPtHisto2D(Form("EtNData%sAntiProton",cutName->Data()),"Number of reconstructed #bar{p}");
 511     CreateEtaPtHisto2D(Form("EtNData%sUnidentified",cutName->Data()),"Number of Reconstructed unidentified particles");
 512
 513     CreateHisto2D(Form("dEdxDataAll%s",cutName->Data()),"dE/dx for all particles","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
 514     CreateHisto2D(Form("dEdxDataPion%s",cutName->Data()),"dE/dx for #pi^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
 515     CreateHisto2D(Form("dEdxDataKaon%s",cutName->Data()),"dE/dx for K^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
 516     CreateHisto2D(Form("dEdxDataProton%s",cutName->Data()),"dE/dx for p(#bar{p})","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
 517     CreateHisto2D(Form("dEdxDataElectron%s",cutName->Data()),"dE/dx for e^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
 518     CreateHisto2D(Form("dEdxDataUnidentified%s",cutName->Data()),"dE/dx for unidentified particles","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
 519   }
 520
 521   Float_t minEt = 0.0;
 522   Float_t maxEt = 100.0;
 523   Int_t nbinsEt = 200;
 524   char histoname[200];
 525   char histotitle[200];
 526   char xtitle[50];
 527   TString *ytitle = new TString("Number of events");
 528   TString *sTPC = new TString("TPC");
 529   TString *sITS = new TString("ITS");
 530   TString *sTPCpt = new TString("0.15");
 531   TString *sITSpt = new TString("0.10");
 532   TString *sPID = new TString("");
 533   TString *sNoPID = new TString("NoPID");
 534   TString *sNoPIDString = new TString(", No PID");
 535   TString *sHadEt = new TString("HadEt");
 536   TString *sRawEt = new TString("RawEt");
 537   TString *sTotEt = new TString("TotEt");
 538   TString *sTotEtString = new TString("total E_{T}");
 539   TString *sHadEtString = new TString("hadronic E_{T}");
 540   TString *sRawEtString = new TString("raw E_{T}");
 541   TString *sFull = new TString("Full");
 542   TString *sEMCAL = new TString("EMCAL");
 543   TString *sPHOS = new TString("PHOS");
 544
 545   for(int tpc = 0;tpc<2;tpc++){
 546     for(int hadet = 0;hadet<3;hadet++){
 547       for(int type = 0;type<3;type++){
 548         for(int pid = 0;pid<2;pid++){
 549           TString *detector;
 550           TString *partid;
 551           TString *et = sHadEt;
 552           TString *acceptance;
 553           TString *ptstring;
 554           TString *partidstring;
 555           TString *etstring = sHadEtString;
 556           if(tpc==1) {detector = sTPC; ptstring = sTPCpt;}
 557           else{detector = sITS; ptstring = sITSpt;}
 558           if(pid==1){partid = sPID; partidstring = sPID;}
 559           else{partid = sNoPID; partidstring = sNoPIDString;}
 560           if(hadet==1) {et = sHadEt; etstring = sHadEtString;}
 561           if(hadet==0){et = sTotEt; etstring = sTotEtString;}
 562           if(hadet==2){et = sRawEt; etstring = sRawEtString;}
 563           switch(type){
 564           case 0:
 565             acceptance = sFull;
 566             break;
 567           case 1:
 568             acceptance = sEMCAL;
 569             break;
 570           case 2:
 571             acceptance = sPHOS;
 572             break;
 573           default:
 574             acceptance = sFull;
 575           }
 576           sprintf(histoname,"Reco%s%sAcceptance%s%s",et->Data(),acceptance->Data(),detector->Data(),partid->Data());
 577           sprintf(histotitle,"Reconstructed %s with %s acceptance for p_{T}>%s GeV/c%s",etstring->Data(),acceptance->Data(),ptstring->Data(),partidstring->Data());
 578           sprintf(xtitle,"Reconstructed %s",etstring->Data());
 579           CreateHisto1D(histoname,histotitle,xtitle,ytitle->Data(),nbinsEt,minEt,maxEt);
 580         }
 581       }
 582     }
 583   }
 584
 585   delete sTPC;
 586   delete sITS;
 587   delete sTPCpt;
 588   delete sITSpt;
 589   delete sPID;
 590   delete sNoPID;
 591   delete sNoPIDString;
 592   delete sHadEt;
 593   delete sTotEt;
 594   delete sTotEtString;
 595   delete sHadEtString;
 596   delete sFull;
 597   delete sEMCAL;
 598   delete sPHOS;
 599
 600 }