PWGLF/totEt/AliAnalysisHadEtMonteCarlo.cxx

   1 //_________________________________________________________________________
   2 //  Utility Class for transverse energy studies, charged hadrons
   3 //  Base class for MC analysis
   4 //  - MC output
   5 // implementation file
   6 //
   7 //Created by Christine Nattrass, Rebecca Scott, Irakli Martashvili
   8 //University of Tennessee at Knoxville
   9 //_________________________________________________________________________
  10 #include "AliAnalysisHadEtMonteCarlo.h"
  11 #include "AliAnalysisEtCuts.h"
  12
  13 #include "AliStack.h"
  14 #include "AliMCEvent.h"
  15 #include "AliESDEvent.h"
  16 #include "AliESDtrackCuts.h"
  17 #include "AliESDpid.h"
  18 #include "AliPID.h"
  19 #include "AliESDtrack.h"
  20 #include "AliVParticle.h"
  21 #include "AliAnalysisTask.h"
  22 #include "AliAnalysisManager.h"
  23 #include "AliAnalysisHadEtReconstructed.h"
  24 #include "AliAnalysisHadEtCorrections.h"
  25 #include "AliAnalysisEtCuts.h"
  26 #include <iostream>
  27 #include "TRandom.h"
  28 #include "AliAnalysisEtCommon.h"
  29 #include "AliCentrality.h"
  30 #include "AliLog.h"
  31 #include "AliPWG0Helper.h"
  32 //class AliPWG0Helper;
  33 //#include "$ALICE_ROOT/PWG0/AliPWG0Helper.h"
  34
  35 using namespace std;
  36
  37 ClassImp(AliAnalysisHadEtMonteCarlo);
  38
  39
  40 Int_t AliAnalysisHadEtMonteCarlo::fgNumSmearWidths = 4;
  41 Float_t AliAnalysisHadEtMonteCarlo::fgSmearWidths[4] = {0.005,0.006,0.007,0.008};
  42
  43 AliAnalysisHadEtMonteCarlo::AliAnalysisHadEtMonteCarlo():AliAnalysisHadEt()
  44                                                         ,fSimPiKPEt(0)
  45                                                         ,fSimHadEt(0)
  46                                                         ,fSimTotEt(0)
  47                                                         ,fSimPiKPEtShouldBeReco(0)
  48                                                         ,fSimPiKPEtShouldBeRecoPi(0)
  49                                                         ,fSimPiKPEtShouldBeRecoK(0)
  50                                                         ,fSimPiKPEtShouldBeRecoP(0)
  51                                                         ,fRunLightweight(0)
  52                                                         ,fInvestigateSmearing(0)
  53                                                         ,fInvestigateFull(0)
  54                                                         ,fInvestigateEMCal(0)
  55                                                         ,fInvestigatePHOS(0)
  56                                                         ,fInvestigatePiKP(0)
  57                                                         ,fRequireITSHits(0)
  58                                                         ,fBaryonEnhancement(0)
  59                                                         ,fUseRecoPt(0)
  60                                                         ,fPtSmearer(0)
  61                                                         ,fHadEtReco(0)
  62 {
  63 }
  64 AliAnalysisHadEtMonteCarlo::~AliAnalysisHadEtMonteCarlo(){//destructor
  65   if(fPtSmearer) delete fPtSmearer;
  66 }
  67
  68 void AliAnalysisHadEtMonteCarlo::ResetEventValues(){//resetting event variables
  69   AliAnalysisHadEt::ResetEventValues();
  70     fSimHadEt=0.0;
  71     fSimTotEt=0.0;
  72     fSimPiKPEt=0.0;
  73 }
  74 Int_t AliAnalysisHadEtMonteCarlo::AnalyseEvent(AliVEvent* ev,AliVEvent* ev2)
  75 { // analyse MC and real event info
  76   FillHisto1D("NEvents",0.5,1);
  77   if(!ev || !ev2){
  78     AliFatal("ERROR: Event does not exist");
  79     return 0;
  80   }
  81   AliMCEvent *mcEvent = dynamic_cast<AliMCEvent*>(ev);
  82   AliESDEvent *realEvent = dynamic_cast<AliESDEvent*>(ev2);
  83   if(!mcEvent || !realEvent){
  84     AliFatal("ERROR: mcEvent or realEvent does not exist");
  85     return 0;
  86   }
  87   AliStack *stack = mcEvent->Stack();
  88   fCentBin= -1;
  89   fGoodEvent = kTRUE;//for p+p collisions if we made it this far we have a good event
  90   if(fDataSet==20100){//If this is Pb+Pb
  91     AliCentrality *centrality = realEvent->GetCentrality();
  92     if(fNCentBins<21) fCentBin= centrality->GetCentralityClass10(fCentralityMethod);
  93     else{ fCentBin= centrality->GetCentralityClass5(fCentralityMethod);}
  94     if(fCentBin ==-1) fGoodEvent = kFALSE;//but for Pb+Pb events we don't want to count events where we did not find a centrality
  95   }
  96   AnalyseEvent(ev);
  97
  98   //for PID
  99   AliESDpid *pID = new AliESDpid();//This is identified as a memory leak in valgrind but I delete this object so I think it may be a problem with AliESDpid.
 100
 101   //=============================================
 102
 103   //Roughly following $ALICE_ROOT/PWG0/dNdEta/AlidNdEtaCorrectionTask
 104
 105   //=============================================TPC&&ITS=============================================
 106   //for investigating momentum smearing
 107   Float_t pTtotalReco = 0.0;
 108   Float_t pTtotalSim = 0.0;
 109   Float_t eTtotalSimAll = 0.0;
 110   Float_t eTtotalReco = 0.0;
 111   Float_t eTtotalRecoEffCorr = 0.0;
 112   Float_t eTtotalRecoEffBkgdCorr = 0.0;
 113   Float_t eTtotalRecoBkgdCorr = 0.0;
 114   Float_t eTtotalRecoUncorr = 0.0;
 115   Float_t eTtotalRecoTotalUncorr = 0.0;
 116   Float_t eTtotalRecoEffCorrPi = 0.0;
 117   Float_t eTtotalRecoEffCorrK = 0.0;
 118   Float_t eTtotalRecoEffCorrP = 0.0;
 119   Float_t eTtotalRecoBkgd = 0.0;
 120   Float_t eTtotalRecoPIDSmeared = 0.0;
 121   Float_t eTtotalAsReconstructed = 0.0;
 122   Float_t eTBkgdAsReconstructed = 0.0;
 123   Float_t eTtotalAsReconstructedPi = 0.0;
 124   Float_t eTtotalAsReconstructedP = 0.0;
 125   Float_t eTtotalAsReconstructedK = 0.0;
 126   Float_t eTtotalSim = 0.0;
 127   Int_t nReco = 0;
 128   TString *strTPC = new TString("TPC");
 129   TString *strITS = new TString("ITS");
 130   TString *strTPCITS = new TString("TPCITS");
 131   Int_t lastcutset = 1;
 132   if(fRequireITSHits) lastcutset = 2;
 133   for(Int_t cutset=0;cutset<=lastcutset;cutset++){
 134     TString *cutName = NULL;
 135     TObjArray* list = NULL;
 136     switch(cutset){
 137     case 0:
 138       cutName = strTPC;
 139       list = fEsdtrackCutsTPC->GetAcceptedTracks(realEvent);
 140       break;
 141     case 1:
 142       cutName = strITS;
 143       list = fEsdtrackCutsITS->GetAcceptedTracks(realEvent);
 144       break;
 145     case 2:
 146       cutName = strTPCITS;
 147       list = fEsdtrackCutsITSTPC->GetAcceptedTracks(realEvent);
 148       break;
 149     default:
 150       cerr<<"Error:  cannot fill histograms!"<<endl;
 151       return -1;
 152     }
 153     Int_t nGoodTracks = list->GetEntries();
 154     for (Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++)
 155       {
 156         AliESDtrack *track = dynamic_cast<AliESDtrack*> (list->At(iTrack));
 157         if (!track)
 158           {
 159             Printf("ERROR: Could not get track %d", iTrack);
 160             continue;
 161           }
 162         else{
 163           Float_t nSigmaPion,nSigmaProton,nSigmaKaon,nSigmaElectron;
 164           pID->MakeTPCPID(track);
 165           pID->MakeITSPID(track);
 166           if(cutset!=1){
 167             nSigmaPion = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kPion));
 168             nSigmaProton = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kProton));
 169             nSigmaKaon = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kKaon));
 170             nSigmaElectron = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kElectron));
 171           }
 172           else{
 173             nSigmaPion = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kPion));
 174             nSigmaProton = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kProton));
 175             nSigmaKaon = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kKaon));
 176             nSigmaElectron = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kElectron));
 177           }
 178 //        bool isPion = (nSigmaPion<3.0 && nSigmaProton>2.0 && nSigmaKaon>2.0);
 179 //        bool isElectron = (nSigmaElectron<2.0 && nSigmaPion>4.0 && nSigmaProton>3.0 && nSigmaKaon>3.0);
 180 //        bool isKaon = (nSigmaPion>3.0 && nSigmaProton>2.0 && nSigmaKaon<2.0);
 181 //        bool isProton = (nSigmaPion>3.0 && nSigmaProton<2.0 && nSigmaKaon>2.0);
 182           bool isPion = (nSigmaPion<3.0 && nSigmaProton>2.0 && nSigmaKaon>2.0);
 183           bool isElectron = (nSigmaElectron<2.0 && nSigmaPion>4.0 && nSigmaProton>3.0 && nSigmaKaon>3.0);
 184           bool isKaon = (nSigmaPion>3.0 && nSigmaProton>3.0 && nSigmaKaon<3.0 && track->Pt()<0.45);
 185           bool isProton = (nSigmaPion>3.0 && nSigmaProton<3.0 && nSigmaKaon>3.0 && track->Pt()<0.9);
 186
 187           bool unidentified = (!isProton && !isKaon && !isElectron && !isPion);
 188           if(cutset==1){//ITS dE/dx identification requires tighter cuts on the tracks and we don't gain much from that so we won't do it
 189             unidentified = true;
 190             isPion=false;
 191             isElectron=false;
 192             isKaon=false;
 193             isProton=false;
 194           }
 195           Float_t dEdx = track->GetTPCsignal();
 196           if(cutset==1) dEdx = track->GetITSsignal();
 197
 198           FillHisto2D(Form("dEdxAll%s",cutName->Data()),track->P(),dEdx,1.0);
 199
 200           UInt_t label = (UInt_t)TMath::Abs(track->GetLabel());
 201           TParticle  *simPart  = stack->Particle(label);
 202           if(!simPart) {
 203             Printf("no MC particle\n");
 204             continue;
 205           }
 206           else{//analysis
 207             if(fInvestigateSmearing && cutset==2){
 208               //calculates what we would measure for the pi/k/p et with background
 209               eTtotalRecoTotalUncorr += Et(simPart);
 210               if(isPion){
 211                 eTtotalRecoEffBkgdCorr += Et(simPart) *fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionPion(track->Pt(),fCentBin) * fHadEtReco->GetCorrections()->GetBackgroundCorrectionTPC(track->Pt());
 212                 eTtotalRecoBkgdCorr += Et(simPart) * fHadEtReco->GetCorrections()->GetBackgroundCorrectionTPC(track->Pt());
 213                 eTtotalRecoTotalUncorr += Et(simPart);
 214               }
 215               if(isProton){
 216                 eTtotalRecoEffBkgdCorr += Et(simPart) *fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionProton(track->Pt(),fCentBin) * fHadEtReco->GetCorrections()->GetBackgroundCorrectionTPC(track->Pt());
 217                 eTtotalRecoBkgdCorr += Et(simPart) * fHadEtReco->GetCorrections()->GetBackgroundCorrectionTPC(track->Pt());
 218               }
 219               if(isKaon){
 220                 eTtotalRecoEffBkgdCorr += Et(simPart) *fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionKaon(track->Pt(),fCentBin) * fHadEtReco->GetCorrections()->GetBackgroundCorrectionTPC(track->Pt());
 221                 eTtotalRecoBkgdCorr += Et(simPart) * fHadEtReco->GetCorrections()->GetBackgroundCorrectionTPC(track->Pt());
 222               }
 223               if(unidentified){
 224                 eTtotalRecoEffBkgdCorr += Et(simPart) *fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionHadron(track->Pt(),fCentBin) * fHadEtReco->GetCorrections()->GetBackgroundCorrectionTPC(track->Pt());
 225                 eTtotalRecoBkgdCorr += Et(simPart) * fHadEtReco->GetCorrections()->GetBackgroundCorrectionTPC(track->Pt());
 226               }
 227               //for calculating et as it's done in the reconstructed data
 228               Float_t corrBkgd=0.0;
 229               Float_t corrNotID=0.0;
 230               Float_t corrNoID=0.0;// = fHadEtReco->GetCorrections()->GetNotIDCorrectionNoPID(track->Pt());
 231               Float_t corrEff = 0.0;
 232               Float_t corrEffNoID = 0.0;
 233               Float_t et = 0.0;
 234               if(cutset==2){//TPC
 235                 corrBkgd = fHadEtReco->GetCorrections()->GetBackgroundCorrectionTPC(track->Pt());
 236                 corrEffNoID = fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionHadron(track->Pt(),fCentBin);
 237                 corrNotID = fHadEtReco->GetCorrections()->GetNotIDConstCorrectionTPC();
 238                 corrNoID = fHadEtReco->GetCorrections()->GetNotIDConstCorrectionTPCNoID();
 239               }
 240               if(cutset==1){//ITS
 241                 corrBkgd = fHadEtReco->GetCorrections()->GetBackgroundCorrectionITS(track->Pt());
 242                 corrEffNoID = fHadEtReco->GetCorrections()->GetITSEfficiencyCorrectionHadron(track->Pt(),fCentBin);
 243                 corrNotID = fHadEtReco->GetCorrections()->GetNotIDConstCorrectionITS();
 244                 corrNoID = fHadEtReco->GetCorrections()->GetNotIDConstCorrectionITSNoID();
 245               }
 246
 247               bool isprimary = stack->IsPhysicalPrimary(label);
 248               if (TMath::Abs(track->Eta()) < fHadEtReco->GetCorrections()->GetEtaCut()){
 249                   if(isPion){
 250                     et = Et(track->P(),track->Theta(),fgPiPlusCode,track->Charge());
 251                     corrEff = fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionPion(track->Pt(),fCentBin);
 252                     if(isprimary){
 253                       eTtotalAsReconstructed += et*corrBkgd*corrEff*corrNotID;
 254                     }
 255                     else{
 256                       eTBkgdAsReconstructed += et*corrBkgd*corrEff*corrNotID;
 257                     }
 258                   }
 259                   if(isKaon){
 260                     et = Et(track->P(),track->Theta(),fgKPlusCode,track->Charge());
 261                     corrEff = fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionKaon(track->Pt(),fCentBin);
 262                     if(isprimary){
 263                       eTtotalAsReconstructed += et*corrBkgd*corrEff*corrNotID;
 264                     }
 265                     else{
 266                       eTBkgdAsReconstructed += et*corrBkgd*corrEff*corrNotID;
 267                     }
 268                   }
 269                   if(isProton){
 270                     et = Et(track->P(),track->Theta(),fgProtonCode,track->Charge());
 271                     corrEff = fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionProton(track->Pt(),fCentBin);
 272                     if(isprimary){
 273                       eTtotalAsReconstructed += et*corrBkgd*corrEff*corrNotID;
 274                     }
 275                     else{
 276                       eTBkgdAsReconstructed += et*corrBkgd*corrEff*corrNotID;
 277                     }
 278                   }
 279                   if(unidentified){
 280                     et = Et(track->P(),track->Theta(),fgPiPlusCode,track->Charge());
 281                     corrEff = fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionHadron(track->Pt(),fCentBin);
 282                     if(isprimary){
 283                       eTtotalAsReconstructed += et*corrBkgd*corrEff*corrNotID;
 284                     }
 285                     else{
 286                       eTBkgdAsReconstructed += et*corrBkgd*corrEff*corrNotID;
 287                     }
 288                   }
 289                   if(!isPion && !isProton && !isKaon && !unidentified){
 290                       eTBkgdAsReconstructed += et*corrBkgd*corrEff*corrNotID;
 291                   }
 292                   Int_t pdgCode =  simPart->GetPDG(0)->PdgCode();
 293                   if(pdgCode==fgPiPlusCode ||pdgCode==fgPiMinusCode){eTtotalAsReconstructedPi+=et*corrBkgd*corrEff*corrNotID;}
 294                   if(pdgCode==fgKPlusCode ||pdgCode==fgKMinusCode){eTtotalAsReconstructedK+=et*corrBkgd*corrEff*corrNotID;}
 295                   if(pdgCode==fgProtonCode ||pdgCode==fgAntiProtonCode){eTtotalAsReconstructedP+=et*corrBkgd*corrEff*corrNotID;}
 296                 }
 297             }
 298
 299             if(cutset==2) eTtotalSimAll += Et(simPart);
 300             if(stack->IsPhysicalPrimary(label)){
 301               if (TMath::Abs(simPart->Eta()) < fHadEtReco->GetCorrections()->GetEtaCut()){
 302                 Int_t pdgCode =  simPart->GetPDG(0)->PdgCode();
 303                 Int_t mypid = 0;
 304                 if(pdgCode==AliAnalysisHadEt::fgPiPlusCode) mypid = 1;
 305                 if(pdgCode==fgProtonCode) mypid = 2;
 306                 if(pdgCode==fgKPlusCode) mypid = 3;
 307                 if(pdgCode==fgEPlusCode) mypid = 4;
 308                 if(pdgCode==fgPiMinusCode) mypid = 1;
 309                 if(pdgCode==fgAntiProtonCode) mypid = 2;
 310                 if(pdgCode==fgKMinusCode) mypid = 3;
 311                 if(pdgCode==fgEMinusCode) mypid = 4;
 312                 bool filled = false;
 313                 //for smearing investigations
 314                 if(fInvestigateSmearing && cutset==2){
 315                   pTtotalReco += simPart->Pt();
 316                   pTtotalSim += track->Pt();
 317                   eTtotalReco += Et(track->P(),track->Theta(),pdgCode,track->Charge());
 318                   eTtotalSim += Et(simPart);
 319                   nReco++;
 320                 }
 321                 //============Charged hadrons===================================
 322                 float myefficiencyCorrEt = 0.0;
 323                 //identified...
 324                 if(isPion){
 325                   if(pdgCode!=fgPiPlusCode && pdgCode!=fgPiMinusCode){
 326                     FillHisto2D(Form("MisidentifiedPIDs%s",cutName->Data()),1,mypid,1);
 327                   }
 328                   float myEt = Et(simPart);
 329                   if(fInvestigateSmearing && cutset==2){
 330                     eTtotalRecoPIDSmeared +=myEt;
 331                     eTtotalRecoEffCorr += myEt *fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionPion(track->Pt(),fCentBin);
 332                     myefficiencyCorrEt = myEt * fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionPion(track->Pt(),fCentBin) ;
 333                   }
 334                   if( !fRunLightweight){
 335                     if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedPiPlus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
 336                     else{ FillHisto2D(Form("EtReconstructed%sIdentifiedPiMinus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
 337                   }
 338                   FillHisto2D(Form("dEdxPion%s",cutName->Data()),track->P(),dEdx,1.0);
 339                 }
 340                 if(isProton){
 341                   if(pdgCode!=fgProtonCode && pdgCode!=fgAntiProtonCode){
 342                     FillHisto2D(Form("MisidentifiedPIDs%s",cutName->Data()),2,mypid,1);
 343                   }
 344                   float myEt = Et(simPart);
 345                   if(fInvestigateSmearing && cutset==2){
 346                     eTtotalRecoPIDSmeared +=myEt;
 347                     eTtotalRecoEffCorr += myEt *fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionProton(track->Pt(),fCentBin);
 348                     myefficiencyCorrEt = myEt * fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionProton(track->Pt(),fCentBin);
 349                   }
 350                   if( !fRunLightweight){
 351                     if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedProton",cutName->Data()),track->Pt(),track->Eta(),myEt);}
 352                     else{ FillHisto2D(Form("EtReconstructed%sIdentifiedAntiProton",cutName->Data()),track->Pt(),track->Eta(),myEt);}
 353                     if(fBaryonEnhancement){
 354                       myEt = myEt*ProtonBaryonEnhancement(track->Pt());
 355                       if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedProtonEnhanced",cutName->Data()),track->Pt(),track->Eta(),myEt);}
 356                       else{ FillHisto2D(Form("EtReconstructed%sIdentifiedAntiProtonEnhanced",cutName->Data()),track->Pt(),track->Eta(),myEt);}
 357                     }
 358                   }
 359                   FillHisto2D(Form("dEdxProton%s",cutName->Data()),track->P(),dEdx,1.0);
 360                 }
 361                 if(isKaon){
 362                   if(pdgCode!=fgKMinusCode && pdgCode!=fgKPlusCode){
 363                     FillHisto2D(Form("MisidentifiedPIDs%s",cutName->Data()),3,mypid,1);
 364                   }
 365                   float myEt = Et(simPart);
 366                   if(fInvestigateSmearing && cutset==2){
 367                     eTtotalRecoPIDSmeared +=myEt;
 368                     eTtotalRecoEffCorr += myEt *fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionKaon(track->Pt(),fCentBin);
 369                     myefficiencyCorrEt = myEt * fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionKaon(track->Pt(),fCentBin);
 370                   }
 371                   if( !fRunLightweight){
 372                     if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedKPlus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
 373                     else{ FillHisto2D(Form("EtReconstructed%sIdentifiedKMinus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
 374                   }
 375                   FillHisto2D(Form("dEdxKaon%s",cutName->Data()),track->P(),dEdx,1.0);
 376                 }
 377                 if(isElectron){
 378                   if(pdgCode!=fgEMinusCode && pdgCode!=fgEPlusCode){
 379                     FillHisto2D(Form("MisidentifiedPIDs%s",cutName->Data()),4,mypid,1);
 380                   }
 381                   if( !fRunLightweight){
 382                     float myEt = Et(simPart);
 383                     if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedEPlus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
 384                     else{ FillHisto2D(Form("EtReconstructed%sIdentifiedEMinus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
 385                   }
 386                   FillHisto2D(Form("dEdxElectron%s",cutName->Data()),track->P(),dEdx,1.0);
 387                 }
 388                 if(unidentified){
 389                   if(pdgCode!=fgEMinusCode && pdgCode!=fgEPlusCode){
 390                     float myEtPi = Et(simPart,fgPionMass);
 391                     float myEtP = Et(simPart,fgProtonMass);
 392                     float myEtK = Et(simPart,fgKaonMass);
 393                     float myEt = Et(simPart);
 394                     if(fInvestigateSmearing && cutset==2){
 395                       eTtotalRecoPIDSmeared +=myEtPi;
 396                       eTtotalRecoEffCorr += myEt *fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionHadron(track->Pt(),fCentBin);
 397                       myefficiencyCorrEt = myEt * fHadEtReco->GetCorrections()->GetTPCEfficiencyCorrectionHadron(track->Pt(),fCentBin);
 398                     }
 399                     if( !fRunLightweight){
 400                       FillHisto2D(Form("EtReconstructed%sUnidentifiedAssumingPion",cutName->Data()),track->Pt(),track->Eta(),myEtPi);
 401                       FillHisto2D(Form("EtReconstructed%sUnidentifiedAssumingProton",cutName->Data()),track->Pt(),track->Eta(),myEtP);
 402                       FillHisto2D(Form("EtReconstructed%sUnidentifiedAssumingKaon",cutName->Data()),track->Pt(),track->Eta(),myEtK);
 403                       FillHisto2D(Form("EtReconstructed%sUnidentified",cutName->Data()),track->Pt(),track->Eta(),myEt);
 404                       FillHisto2D(Form("EtNReconstructed%sUnidentified",cutName->Data()),track->Pt(),track->Eta(),1.0);
 405                       if(pdgCode == fgPiPlusCode||pdgCode == fgPiMinusCode){
 406                         FillHisto2D(Form("EtReconstructed%sUnidentifiedPionAssumingPion",cutName->Data()),track->Pt(),track->Eta(),myEtPi);
 407                         FillHisto2D(Form("EtReconstructed%sUnidentifiedPionAssumingProton",cutName->Data()),track->Pt(),track->Eta(),myEtP);
 408                         FillHisto2D(Form("EtReconstructed%sUnidentifiedPionAssumingKaon",cutName->Data()),track->Pt(),track->Eta(),myEtK);
 409                         FillHisto2D(Form("EtReconstructed%sUnidentifiedPion",cutName->Data()),track->Pt(),track->Eta(),myEt);
 410                         FillHisto2D(Form("EtNReconstructed%sUnidentifiedPion",cutName->Data()),track->Pt(),track->Eta(),1.0);
 411                       }
 412                       if(pdgCode == fgKPlusCode||pdgCode == fgKMinusCode){
 413                         FillHisto2D(Form("EtReconstructed%sUnidentifiedKaonAssumingPion",cutName->Data()),track->Pt(),track->Eta(),myEtPi);
 414                         FillHisto2D(Form("EtReconstructed%sUnidentifiedKaonAssumingProton",cutName->Data()),track->Pt(),track->Eta(),myEtP);
 415                         FillHisto2D(Form("EtReconstructed%sUnidentifiedKaonAssumingKaon",cutName->Data()),track->Pt(),track->Eta(),myEtK);
 416                         FillHisto2D(Form("EtReconstructed%sUnidentifiedKaon",cutName->Data()),track->Pt(),track->Eta(),myEt);
 417                         FillHisto2D(Form("EtNReconstructed%sUnidentifiedKaon",cutName->Data()),track->Pt(),track->Eta(),1.0);
 418                       }
 419                       if(pdgCode == fgProtonCode||pdgCode == fgAntiProtonCode){
 420                         FillHisto2D(Form("EtReconstructed%sUnidentifiedProtonAssumingPion",cutName->Data()),track->Pt(),track->Eta(),myEtPi);
 421                         FillHisto2D(Form("EtReconstructed%sUnidentifiedProtonAssumingProton",cutName->Data()),track->Pt(),track->Eta(),myEtP);
 422                         FillHisto2D(Form("EtReconstructed%sUnidentifiedProtonAssumingKaon",cutName->Data()),track->Pt(),track->Eta(),myEtK);
 423                         FillHisto2D(Form("EtReconstructed%sUnidentifiedProton",cutName->Data()),track->Pt(),track->Eta(),myEt);
 424                         FillHisto2D(Form("EtNReconstructed%sUnidentifiedProton",cutName->Data()),track->Pt(),track->Eta(),1.0);
 425                         if(fBaryonEnhancement){
 426                           myEt = myEt*ProtonBaryonEnhancement(track->Pt());
 427                           FillHisto2D(Form("EtReconstructed%sUnidentifiedProtonAssumingPionEnhanced",cutName->Data()),track->Pt(),track->Eta(),myEtPi);
 428                           FillHisto2D(Form("EtReconstructed%sUnidentifiedProtonEnhanced",cutName->Data()),track->Pt(),track->Eta(),myEt);
 429                           FillHisto2D(Form("EtNReconstructed%sUnidentifiedProtonEnhanced",cutName->Data()),track->Pt(),track->Eta(),1.0);
 430                         }
 431                       }
 432                     }
 433                   }
 434                   FillHisto2D(Form("dEdxUnidentified%s",cutName->Data()),track->P(),dEdx,1.0);
 435                   FillHisto1D(Form("UnidentifiedPIDs%s",cutName->Data()),mypid,1);
 436                 }
 437                 //...simulated
 438                 float myEtSim = Et(simPart);
 439                 float myEtReco = 0.0;
 440                 if(pdgCode == fgPiPlusCode){
 441                   float myEt = Et(simPart);
 442                   float myEtP = Et(simPart,fgProtonMass);
 443                   float myEtK = Et(simPart,fgKaonMass);
 444                   myEtReco = Et(track->P(),track->Theta(),fgPiPlusCode,track->Charge());
 445                   float pT = simPart->Pt();
 446                   float eta = simPart->Eta();
 447                   if(fInvestigateSmearing && cutset==2){
 448                     eTtotalRecoEffCorrPi+=myefficiencyCorrEt;
 449                     eTtotalRecoUncorr +=myEt;
 450                   }
 451                   if(fUseRecoPt){//Then we switch the pT and the Et
 452                     myEt = myEtReco;
 453                     pT = track->Pt();
 454                     eta = track->Eta();
 455                   }
 456                   if( !fRunLightweight){
 457                     FillHisto2D(Form("EtReconstructed%sPiPlus",cutName->Data()),pT,eta,myEt);
 458                     FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 459                     FillHisto2D(Form("EtNReconstructed%sPiPlus",cutName->Data()),pT,eta,myEt);
 460                     FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 461                     if(fCentBin>=0){//if a centrality bin was defined
 462                       FillHisto2D(Form("EtNReconstructed%sPiPlusCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 463                       FillHisto2D(Form("EtNReconstructed%sChargedHadronCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 464                     }
 465                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),pT,eta,myEt);
 466                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingProton",cutName->Data()),pT,eta,myEtP);
 467                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingKaon",cutName->Data()),pT,eta,myEtK);
 468                     FillHisto2D(Form("EtReconstructed%sPiPlusAssumingKaon",cutName->Data()),pT,eta,myEtK);
 469                     FillHisto2D(Form("EtReconstructed%sPiPlusAssumingProton",cutName->Data()),pT,eta,myEtP);
 470                   }
 471                   filled = true;
 472                 }
 473                 if(pdgCode == fgPiMinusCode){
 474                   float myEt = Et(simPart);
 475                   float myEtP = Et(simPart,fgProtonMass);
 476                   float myEtK = Et(simPart,fgKaonMass);
 477                   myEtReco = Et(track->P(),track->Theta(),fgPiMinusCode,track->Charge());
 478                   float pT = simPart->Pt();
 479                   float eta = simPart->Eta();
 480                   if(fInvestigateSmearing && cutset==2){
 481                     eTtotalRecoEffCorrPi+=myefficiencyCorrEt;
 482                     eTtotalRecoUncorr +=myEt;
 483                   }
 484                   if(fUseRecoPt){//Then we switch the pT and the Et
 485                     myEt = myEtReco;
 486                     pT = track->Pt();
 487                     eta = track->Eta();
 488                   }
 489                   if( !fRunLightweight){
 490                     FillHisto2D(Form("EtReconstructed%sPiMinus",cutName->Data()),pT,eta,myEt);
 491                     FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 492                     FillHisto2D(Form("EtNReconstructed%sPiMinus",cutName->Data()),pT,eta,myEt);
 493                     FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 494                     if(fCentBin>=0){//if a centrality bin was defined
 495                       FillHisto2D(Form("EtNReconstructed%sPiMinusCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 496                       FillHisto2D(Form("EtNReconstructed%sChargedHadronCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 497                     }
 498                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),pT,eta,myEt);
 499                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingProton",cutName->Data()),pT,eta,myEtP);
 500                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingKaon",cutName->Data()),pT,eta,myEtK);
 501                     FillHisto2D(Form("EtReconstructed%sPiMinusAssumingKaon",cutName->Data()),pT,eta,myEtK);
 502                     FillHisto2D(Form("EtReconstructed%sPiMinusAssumingProton",cutName->Data()),pT,eta,myEtP);
 503                   }
 504                   filled = true;
 505                 }
 506                 if(pdgCode == fgKPlusCode){
 507                   float myEt = Et(simPart);
 508                   float myEtPi = Et(simPart,fgPionMass);
 509                   float myEtP = Et(simPart,fgProtonMass);
 510                   myEtReco = Et(track->P(),track->Theta(),fgKPlusCode,track->Charge());
 511                   float pT = simPart->Pt();
 512                   float eta = simPart->Eta();
 513                   if(fUseRecoPt){//Then we switch the pT and the Et
 514                     myEt = myEtReco;
 515                     pT = track->Pt();
 516                     eta = track->Eta();
 517                   }
 518                   if(fInvestigateSmearing && cutset==2){
 519                     eTtotalRecoEffCorrK+=myefficiencyCorrEt;
 520                     eTtotalRecoUncorr +=myEt;
 521                   }
 522                   if( !fRunLightweight){
 523                     FillHisto2D(Form("EtReconstructed%sKPlus",cutName->Data()),pT,eta,myEt);
 524                     FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 525                     FillHisto2D(Form("EtNReconstructed%sKPlus",cutName->Data()),pT,eta,myEt);
 526                     FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 527                     if(fCentBin>=0){//if a centrality bin was defined
 528                       FillHisto2D(Form("EtNReconstructed%sKPlusCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 529                       FillHisto2D(Form("EtNReconstructed%sChargedHadronCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 530                     }
 531                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),pT,eta,myEtPi);
 532                     FillHisto2D(Form("EtReconstructed%sKPlusAssumingPion",cutName->Data()),pT,eta,myEtPi);
 533                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingKaon",cutName->Data()),pT,eta,myEt);
 534                     FillHisto2D(Form("EtReconstructed%sKPlusAssumingKaon",cutName->Data()),pT,eta,myEt);
 535                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingProton",cutName->Data()),pT,eta,myEtP);
 536                     FillHisto2D(Form("EtReconstructed%sKPlusAssumingProton",cutName->Data()),pT,eta,myEtP);
 537                   }
 538                   filled = true;
 539                 }
 540                 if(pdgCode == fgKMinusCode){
 541                   float myEt = Et(simPart);
 542                   float myEtPi = Et(simPart,fgPionMass);
 543                   float myEtP = Et(simPart,fgProtonMass);
 544                   myEtReco = Et(track->P(),track->Theta(),fgKMinusCode,track->Charge());
 545                   float pT = simPart->Pt();
 546                   float eta = simPart->Eta();
 547                   if(fUseRecoPt){//Then we switch the pT and the Et
 548                     myEt = myEtReco;
 549                     pT = track->Pt();
 550                     eta = track->Eta();
 551                   }
 552                   if(fInvestigateSmearing && cutset==2){
 553                     eTtotalRecoEffCorrK+=myefficiencyCorrEt;
 554                     eTtotalRecoUncorr +=myEt;
 555                   }
 556                   if( !fRunLightweight){
 557                     FillHisto2D(Form("EtReconstructed%sKMinus",cutName->Data()),pT,eta,myEt);
 558                     FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 559                     FillHisto2D(Form("EtNReconstructed%sKMinus",cutName->Data()),pT,eta,myEt);
 560                     FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 561                     if(fCentBin>=0){//if a centrality bin was defined
 562                       FillHisto2D(Form("EtNReconstructed%sKMinusCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 563                       FillHisto2D(Form("EtNReconstructed%sChargedHadronCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 564                     }
 565                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),pT,eta,myEtPi);
 566                     FillHisto2D(Form("EtReconstructed%sKMinusAssumingPion",cutName->Data()),pT,eta,myEtPi);
 567                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingKaon",cutName->Data()),pT,eta,myEt);
 568                     FillHisto2D(Form("EtReconstructed%sKMinusAssumingKaon",cutName->Data()),pT,eta,myEt);
 569                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingProton",cutName->Data()),pT,eta,myEtP);
 570                     FillHisto2D(Form("EtReconstructed%sKMinusAssumingProton",cutName->Data()),pT,eta,myEtP);
 571                   }
 572                   filled = true;
 573                 }
 574                 if(pdgCode == fgProtonCode){
 575                   float myEt = Et(simPart);
 576                   float myEtPi = Et(simPart,fgPionMass);
 577                   float myEtK = Et(simPart,fgKaonMass);
 578                   myEtReco = Et(track->P(),track->Theta(),fgProtonCode,track->Charge());
 579                   float pT = simPart->Pt();
 580                   float eta = simPart->Eta();
 581                   if(fUseRecoPt){//Then we switch the pT and the Et
 582                     myEt = myEtReco;
 583                     pT = track->Pt();
 584                     eta = track->Eta();
 585                   }
 586                   if(fInvestigateSmearing && cutset==2){
 587                     eTtotalRecoEffCorrP+=myefficiencyCorrEt;
 588                     eTtotalRecoUncorr +=myEt;
 589                   }
 590                   if( !fRunLightweight){
 591                     FillHisto2D(Form("EtReconstructed%sProton",cutName->Data()),pT,eta,myEt);
 592                     FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 593                     FillHisto2D(Form("EtNReconstructed%sProton",cutName->Data()),pT,eta,myEt);
 594                     FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 595                     if(fCentBin>=0){//if a centrality bin was defined
 596                       FillHisto2D(Form("EtNReconstructed%sProtonCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 597                       FillHisto2D(Form("EtNReconstructed%sChargedHadronCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 598                     }
 599                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),pT,eta,myEtPi);
 600                     FillHisto2D(Form("EtReconstructed%sProtonAssumingPion",cutName->Data()),pT,eta,myEtPi);
 601                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingKaon",cutName->Data()),pT,eta,myEtK);
 602                     FillHisto2D(Form("EtReconstructed%sProtonAssumingKaon",cutName->Data()),pT,eta,myEtK);
 603                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingProton",cutName->Data()),pT,eta,myEt);
 604                     FillHisto2D(Form("EtReconstructed%sProtonAssumingProton",cutName->Data()),pT,eta,myEt);
 605                   }
 606                   filled = true;
 607
 608                   if( !fRunLightweight){
 609                     if(fBaryonEnhancement){
 610                       float enhancement = ProtonBaryonEnhancement(track->Pt());
 611                       FillHisto2D(Form("EtReconstructed%sProtonEnhanced",cutName->Data()),pT,eta,myEt*enhancement);
 612                       FillHisto2D(Form("EtNReconstructed%sProtonEnhanced",cutName->Data()),pT,eta,myEt*enhancement);
 613                       FillHisto2D(Form("EtReconstructed%sProtonAssumingPionEnhanced",cutName->Data()),pT,eta,myEtPi*enhancement);
 614                     }
 615                   }
 616
 617                 }
 618                 if(pdgCode == fgAntiProtonCode){
 619                   float myEt = Et(simPart);
 620                   float myEtPi = Et(simPart,fgPionMass);
 621                   float myEtK = Et(simPart,fgKaonMass);
 622                   myEtReco = Et(track->P(),track->Theta(),fgAntiProtonCode,track->Charge());
 623                   float pT = simPart->Pt();
 624                   float eta = simPart->Eta();
 625                   if(fUseRecoPt){//Then we switch the pT and the Et
 626                     myEt = myEtReco;
 627                     pT = track->Pt();
 628                     eta = track->Eta();
 629                   }
 630                   if(fInvestigateSmearing && cutset==2){
 631                     eTtotalRecoEffCorrP+=myefficiencyCorrEt;
 632                     eTtotalRecoUncorr +=myEt;
 633                   }
 634                   if( !fRunLightweight){
 635                     FillHisto2D(Form("EtReconstructed%sAntiProton",cutName->Data()),pT,eta,myEt);
 636                     FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 637                     FillHisto2D(Form("EtNReconstructed%sAntiProton",cutName->Data()),pT,eta,myEt);
 638                     FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),pT,eta,myEt);
 639                     if(fCentBin>=0){//if a centrality bin was defined
 640                       FillHisto2D(Form("EtNReconstructed%sAntiProtonCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 641                       FillHisto2D(Form("EtNReconstructed%sChargedHadronCB%i",cutName->Data(),fCentBin),pT,eta,myEt);
 642                     }
 643                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),pT,eta,myEtPi);
 644                     FillHisto2D(Form("EtReconstructed%sAntiProtonAssumingPion",cutName->Data()),pT,eta,myEtPi);
 645                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingKaon",cutName->Data()),pT,eta,myEtK);
 646                     FillHisto2D(Form("EtReconstructed%sAntiProtonAssumingKaon",cutName->Data()),pT,eta,myEtK);
 647                     FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingProton",cutName->Data()),pT,eta,myEt);
 648                     FillHisto2D(Form("EtReconstructed%sAntiProtonAssumingProton",cutName->Data()),pT,eta,myEt);
 649                   }
 650                   filled = true;
 651                   if( !fRunLightweight){
 652                     if(fBaryonEnhancement){
 653                       float enhancement = ProtonBaryonEnhancement(track->Pt());
 654                       FillHisto2D(Form("EtReconstructed%sAntiProtonEnhanced",cutName->Data()),pT,eta,myEt*enhancement);
 655                       FillHisto2D(Form("EtNReconstructed%sAntiProtonEnhanced",cutName->Data()),pT,eta,myEt*enhancement);
 656                       FillHisto2D(Form("EtReconstructed%sAntiProtonAssumingPionEnhanced",cutName->Data()),pT,eta,myEtPi*enhancement);
 657                     }
 658                   }
 659                 }
 660                 if(pdgCode == fgEPlusCode){
 661                   if( !fRunLightweight){
 662                     float myEt = Et(simPart);
 663                     FillHisto2D(Form("EtReconstructed%sEPlus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
 664                     if(!isElectron || unidentified){
 665                       float myEtPi = Et(simPart,fgPionMass);
 666                       FillHisto2D(Form("EtReconstructed%sMisidentifiedElectrons",cutName->Data()),track->Pt(),track->Eta(),myEtPi);
 667                     }
 668                   }
 669                   filled = true;
 670                 }
 671                 if(pdgCode == fgEMinusCode){
 672                   if( !fRunLightweight){
 673                     if(!isElectron || unidentified){
 674                       float myEtPi = Et(simPart,fgPionMass);
 675                       FillHisto2D(Form("EtReconstructed%sMisidentifiedElectrons",cutName->Data()),track->Pt(),track->Eta(),myEtPi);
 676                     }
 677                     float myEt = Et(simPart);
 678                     FillHisto2D(Form("EtReconstructed%sEMinus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
 679                   }
 680                   filled = true;
 681                 }
 682                 if( !fRunLightweight){
 683                   if(myEtReco>0.0){FillHisto2D(Form("ETresolution%s",cutName->Data()),myEtReco,(myEtSim-myEtReco)/myEtReco,1.0);}
 684                   if(track->Pt()>0.0){FillHisto2D(Form("pTresolution%s",cutName->Data()),track->Pt(),(simPart->Pt() - track->Pt())/track->Pt(),1.0);}
 685                   if(track->P()>0.0){FillHisto2D(Form("presolution%s",cutName->Data()),track->P(),(simPart->P() - track->P())/track->P(),1.0);}
 686                   FillHisto1D(Form("pTsim%s",cutName->Data()),simPart->Pt(),1.0);
 687                   FillHisto1D(Form("pTrec%s",cutName->Data()),track->Pt(),1.0);
 688                   if(fCentBin!=-1){
 689                     FillHisto1D(Form("pTsim%sCB%i",cutName->Data(),fCentBin),simPart->Pt(),1.0);
 690                     FillHisto1D(Form("pTrec%sCB%i",cutName->Data(),fCentBin),track->Pt(),1.0);
 691                   }
 692                 }
 693               }
 694             }
 695             else{//not a primary - we're after V0 daughters!
 696               bool written = false;
 697               //now, what is the et we would measure for this?  Since this is the relevant et.
 698               float myrecoEt = 0;
 699               if(isPion || unidentified) myrecoEt = Et(track->P(),track->Theta(),fgPiPlusCode,track->Charge());
 700               if(isProton) myrecoEt = Et(track->P(),track->Theta(),fgProtonCode,track->Charge());
 701               if(isKaon) myrecoEt = Et(track->P(),track->Theta(),fgKPlusCode,track->Charge());
 702               if (TMath::Abs(simPart->Eta()) < fHadEtReco->GetCorrections()->GetEtaCut()){
 703                 TParticle *mom = stack->Particle(simPart->GetFirstMother());
 704                 if(mom){
 705                   TParticlePDG *pc = mom->GetPDG(0);
 706                   if(pc){
 707                     Int_t pdgCode =  mom->GetPDG(0)->PdgCode();
 708                     if(pdgCode == fgLambdaCode){
 709                       written = true;
 710                       float myEt = Et(simPart);
 711                       float pT = simPart->Pt();
 712                       float eta = simPart->Eta();
 713                       eTtotalRecoBkgd+=myEt;
 714                       if(fUseRecoPt){//Then we switch the pT and the Et
 715                         myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 716                         pT = track->Pt();
 717                         eta = track->Eta();
 718                       }
 719                       if( !fRunLightweight){
 720                         FillHisto2D(Form("EtReconstructed%sLambdaDaughters",cutName->Data()),pT,eta,myrecoEt);
 721                         Float_t weight = LambdaWeight(mom->Pt());
 722                         if(fBaryonEnhancement){
 723                           float enhancement = ProtonBaryonEnhancement(track->Pt());
 724                           weight = weight*enhancement;
 725                         }
 726                         FillHisto2D(Form("EtReconstructed%sLambdaDaughtersReweighted",cutName->Data()),pT,eta,myrecoEt*weight);
 727                       }
 728                     }
 729                     if(pdgCode == fgAntiLambdaCode){
 730                       written = true;
 731                       float myEt = Et(simPart);
 732                       float pT = simPart->Pt();
 733                       float eta = simPart->Eta();
 734                       eTtotalRecoBkgd+=myEt;
 735                       if(fUseRecoPt){//Then we switch the pT and the Et
 736                         myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 737                         pT = track->Pt();
 738                         eta = track->Eta();
 739                       }
 740                       if( !fRunLightweight){
 741                         FillHisto2D(Form("EtReconstructed%sAntiLambdaDaughters",cutName->Data()),pT,eta,myrecoEt);
 742                         Float_t weight = AntiLambdaWeight(mom->Pt());
 743                         if(fBaryonEnhancement){
 744                           float enhancement = ProtonBaryonEnhancement(track->Pt());
 745                           weight = weight*enhancement;
 746                         }
 747                         FillHisto2D(Form("EtReconstructed%sAntiLambdaDaughtersReweighted",cutName->Data()),pT,eta,myrecoEt*weight);
 748                       }
 749                     }
 750                     if(pdgCode == fgK0SCode || pdgCode == fgK0LCode || pdgCode == fgKPlusCode || pdgCode == fgKMinusCode){//actually get all kaon daughters
 751                       written = true;
 752                       float myEt = Et(simPart);
 753                       float pT = simPart->Pt();
 754                       float eta = simPart->Eta();
 755                       eTtotalRecoBkgd+=myEt;
 756                       if(fUseRecoPt){//Then we switch the pT and the Et
 757                         myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 758                         pT = track->Pt();
 759                         eta = track->Eta();
 760                       }
 761                       if( !fRunLightweight){
 762                         FillHisto2D(Form("EtReconstructed%sK0SDaughters",cutName->Data()),pT,eta,myEt);
 763                         Float_t weight = K0Weight(mom->Pt());
 764                         FillHisto2D(Form("EtReconstructed%sK0SDaughtersReweighted",cutName->Data()),pT,eta,myrecoEt*weight);
 765                       }
 766                     }
 767                     if(pdgCode == fgXiCode){
 768                       written = true;
 769                       float myEt = Et(simPart);
 770                       float pT = simPart->Pt();
 771                       float eta = simPart->Eta();
 772                       eTtotalRecoBkgd+=myEt;
 773                       if(fUseRecoPt){//Then we switch the pT and the Et
 774                         myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 775                         pT = track->Pt();
 776                         eta = track->Eta();
 777                       }
 778                       if( !fRunLightweight){
 779                         FillHisto2D(Form("EtReconstructed%sXiDaughters",cutName->Data()),pT,eta,myrecoEt);
 780                       }
 781                     }
 782                     if(pdgCode == fgAntiXiCode){
 783                       written = true;
 784                       float myEt = Et(simPart);
 785                       float pT = simPart->Pt();
 786                       float eta = simPart->Eta();
 787                       eTtotalRecoBkgd+=myEt;
 788                       if(fUseRecoPt){//Then we switch the pT and the Et
 789                         myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 790                         pT = track->Pt();
 791                         eta = track->Eta();
 792                       }
 793                       if( !fRunLightweight){
 794                         FillHisto2D(Form("EtReconstructed%sAntiXiDaughters",cutName->Data()),pT,eta,myrecoEt);
 795                       }
 796                     }
 797                     if(pdgCode == fgOmegaCode){
 798                       written = true;
 799                       float myEt = Et(simPart);
 800                       float pT = simPart->Pt();
 801                       float eta = simPart->Eta();
 802                       eTtotalRecoBkgd+=myEt;
 803                       if(fUseRecoPt){//Then we switch the pT and the Et
 804                         myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 805                         pT = track->Pt();
 806                         eta = track->Eta();
 807                       }
 808                       if( !fRunLightweight){
 809                         FillHisto2D(Form("EtReconstructed%sOmegaDaughters",cutName->Data()),pT,eta,myrecoEt);
 810                       }
 811                     }
 812                     if(pdgCode == fgXiCode){
 813                       written = true;
 814                       float myEt = Et(simPart);
 815                       float pT = simPart->Pt();
 816                       float eta = simPart->Eta();
 817                       eTtotalRecoBkgd+=myEt;
 818                       if(fUseRecoPt){//Then we switch the pT and the Et
 819                         myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 820                         pT = track->Pt();
 821                         eta = track->Eta();
 822                       }
 823                       if( !fRunLightweight){
 824                         FillHisto2D(Form("EtReconstructed%sAntiOmegaDaughters",cutName->Data()),pT,eta,myrecoEt);
 825                       }
 826                     }
 827
 828                     if(mom->GetFirstMother()>0){
 829                       TParticle *grandma = stack->Particle(mom->GetFirstMother());
 830                       if(grandma){
 831                         Int_t pdgCodeMom =  mom->GetPDG(0)->PdgCode();
 832                         if(pdgCodeMom==fgPiPlusCode || pdgCodeMom==fgPiMinusCode || pdgCodeMom==fgProtonCode ||pdgCodeMom==fgAntiProtonCode || pdgCodeMom==fgKPlusCode || pdgCode==fgKMinusCode){
 833                           Int_t pdgCodeGrandma =  grandma->GetPDG(0)->PdgCode();
 834
 835                           if(pdgCodeGrandma == fgXiCode){
 836                             written = true;
 837                             float myEt = Et(simPart);
 838                             float pT = simPart->Pt();
 839                             float eta = simPart->Eta();
 840                             eTtotalRecoBkgd+=myEt;
 841                             if(fUseRecoPt){//Then we switch the pT and the Et
 842                               myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 843                               pT = track->Pt();
 844                               eta = track->Eta();
 845                             }
 846                             if( !fRunLightweight){
 847                               FillHisto2D(Form("EtReconstructed%sXiDaughters",cutName->Data()),pT,eta,myrecoEt);
 848                             }
 849                           }
 850                           if(pdgCodeGrandma == fgAntiXiCode){
 851                             written = true;
 852                             float myEt = Et(simPart);
 853                             float pT = simPart->Pt();
 854                             float eta = simPart->Eta();
 855                             eTtotalRecoBkgd+=myEt;
 856                             if(fUseRecoPt){//Then we switch the pT and the Et
 857                               myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 858                               pT = track->Pt();
 859                               eta = track->Eta();
 860                             }
 861                             if( !fRunLightweight){
 862                               FillHisto2D(Form("EtReconstructed%sAntiXiDaughters",cutName->Data()),pT,eta,myrecoEt);
 863                             }
 864                           }
 865                           if(pdgCodeGrandma == fgOmegaCode){
 866                             written = true;
 867                             float myEt = Et(simPart);
 868                             float pT = simPart->Pt();
 869                             float eta = simPart->Eta();
 870                             eTtotalRecoBkgd+=myEt;
 871                             if(fUseRecoPt){//Then we switch the pT and the Et
 872                               myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 873                               pT = track->Pt();
 874                               eta = track->Eta();
 875                             }
 876                             if( !fRunLightweight){
 877                               FillHisto2D(Form("EtReconstructed%sOmegaDaughters",cutName->Data()),pT,eta,myrecoEt);
 878                             }
 879                           }
 880                           if(pdgCodeGrandma == fgXiCode){
 881                             written = true;
 882                             float myEt = Et(simPart);
 883                             float pT = simPart->Pt();
 884                             float eta = simPart->Eta();
 885                             eTtotalRecoBkgd+=myEt;
 886                             if(fUseRecoPt){//Then we switch the pT and the Et
 887                               myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 888                               pT = track->Pt();
 889                               eta = track->Eta();
 890                             }
 891                             if( !fRunLightweight){
 892                               FillHisto2D(Form("EtReconstructed%sAntiOmegaDaughters",cutName->Data()),pT,eta,myrecoEt);
 893                             }
 894                           }
 895
 896                         }
 897                       }
 898                     }
 899                     if(!written){
 900                       int mycode = simPart->GetPDG(0)->PdgCode();
 901                       if( (pdgCode == fgGammaCode || pdgCode == fgPi0Code) && (mycode==fgEPlusCode||mycode==fgEMinusCode)){
 902                         written = true;
 903                         float myEt = Et(simPart);
 904                         float pT = simPart->Pt();
 905                         float eta = simPart->Eta();
 906                         eTtotalRecoBkgd+=myEt;
 907                         if(fUseRecoPt){//Then we switch the pT and the Et
 908                           myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 909                           pT = track->Pt();
 910                           eta = track->Eta();
 911                         }
 912                         if( !fRunLightweight){
 913                           FillHisto2D(Form("EtReconstructed%sConversionElectrons",cutName->Data()),pT,eta,myrecoEt);
 914                         }
 915                       }
 916                       if(mycode==fgMuPlusCode || mycode==fgMuMinusCode){
 917                         written = true;
 918                         float myEt = Et(simPart);
 919                         float pT = simPart->Pt();
 920                         float eta = simPart->Eta();
 921                         eTtotalRecoBkgd+=myEt;
 922                         if(fUseRecoPt){//Then we switch the pT and the Et
 923                           myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 924                           pT = track->Pt();
 925                           eta = track->Eta();
 926                         }
 927                         if( !fRunLightweight){
 928                           FillHisto2D(Form("EtReconstructed%sSecondaryMuons",cutName->Data()),pT,eta,myrecoEt);
 929                         }
 930                       }
 931                       if(mycode==fgPiPlusCode || mycode==fgPiMinusCode){
 932                         written = true;
 933                         float myEt = Et(simPart);
 934                         float pT = simPart->Pt();
 935                         float eta = simPart->Eta();
 936                         eTtotalRecoBkgd+=myEt;
 937                         if(fUseRecoPt){//Then we switch the pT and the Et
 938                           myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 939                           pT = track->Pt();
 940                           eta = track->Eta();
 941                         }
 942                         if( !fRunLightweight){
 943                           FillHisto2D(Form("EtReconstructed%sSecondaryPions",cutName->Data()),pT,eta,myrecoEt);
 944                         }
 945                       }
 946                       if(mycode==fgAntiProtonCode || mycode==fgProtonCode){
 947                         written = true;
 948                         float myEt = Et(simPart);
 949                         float pT = simPart->Pt();
 950                         float eta = simPart->Eta();
 951                         eTtotalRecoBkgd+=myEt;
 952                         if(fUseRecoPt){//Then we switch the pT and the Et
 953                           myEt = Et(track->P(),track->Theta(),simPart->GetPDG(0)->PdgCode(),track->Charge());
 954                           pT = track->Pt();
 955                           eta = track->Eta();
 956                         }
 957                         if( !fRunLightweight){
 958                           FillHisto2D(Form("EtReconstructed%sSecondaryProtons",cutName->Data()),pT,eta,myrecoEt);
 959                         }
 960                       }
 961                       //if(!written) cout<<"I was not counted in the background and I am a "<<simPart->GetName()<<" and my mother is a "<<mom->GetName()<<endl;
 962                     }
 963                   }
 964                   else{cout<<"No particle code!! 657"<<endl;}
 965                 }
 966                 else{cout<<"No mother particle!! 658"<<endl;}
 967               }
 968             }
 969           }
 970
 971         }
 972       }
 973     delete list;
 974   }
 975   if(fInvestigateSmearing && !fRunLightweight){
 976     if(fSimPiKPEtShouldBeReco>0.0) FillHisto2D("SimPiKPEtMinusSimEffCorrRecoOnly",fSimPiKPEtShouldBeReco,(fSimPiKPEtShouldBeReco-eTtotalRecoEffCorr)/fSimPiKPEtShouldBeReco,1.0);
 977     if(fSimPiKPEtShouldBeReco>0.0) FillHisto2D("SimPiKPEtMinusSimEffBkgdCorrRecoOnly",fSimPiKPEtShouldBeReco,(fSimPiKPEtShouldBeReco-eTtotalRecoEffBkgdCorr)/fSimPiKPEtShouldBeReco,1.0);
 978     if(fSimPiKPEtShouldBeRecoPi>0.0) FillHisto2D("SimPiKPEtMinusSimEffCorrRecoPiOnly",fSimPiKPEtShouldBeRecoPi,(fSimPiKPEtShouldBeRecoPi-eTtotalRecoEffCorrPi)/fSimPiKPEtShouldBeRecoPi,1.0);
 979     if(fSimPiKPEtShouldBeRecoP>0.0) FillHisto2D("SimPiKPEtMinusSimEffCorrRecoPOnly",fSimPiKPEtShouldBeRecoP,(fSimPiKPEtShouldBeRecoP-eTtotalRecoEffCorrP)/fSimPiKPEtShouldBeRecoP,1.0);
 980     if(fSimPiKPEtShouldBeRecoK>0.0) FillHisto2D("SimPiKPEtMinusSimEffCorrRecoKOnly",fSimPiKPEtShouldBeRecoK,(fSimPiKPEtShouldBeRecoK-eTtotalRecoEffCorrK)/fSimPiKPEtShouldBeRecoK,1.0);
 981     if(eTtotalSim>0.0) FillHisto2D("SimPiKPEtMinusSimAllCorrSmearedRecoOnly",eTtotalSim,(eTtotalSim-eTtotalSimAll+eTtotalRecoBkgd)/eTtotalSim,1.0);
 982     if(eTtotalRecoTotalUncorr>0.0) FillHisto2D("SimPiKPEtMeasMinusEtRealPiKP",eTtotalRecoTotalUncorr,(eTtotalRecoTotalUncorr-eTtotalRecoUncorr)/eTtotalRecoTotalUncorr,1.0);
 983     if(eTtotalSim>0.0) FillHisto2D("SimPiKPEtMinusSimAllSmearedRecoOnly",eTtotalSim,(eTtotalSim-eTtotalSimAll)/eTtotalSim,1.0);
 984     if(eTtotalSim>0.0) FillHisto2D("SimPiKPEtMinusSimPIDSmearedRecoOnly",eTtotalSim,(eTtotalSim-eTtotalRecoPIDSmeared*1.01)/eTtotalSim,1.0);
 985     if(eTtotalSim>0.0) FillHisto2D("SimPiKPEtMinusSimSmearedRecoOnly",eTtotalSim,(eTtotalSim-eTtotalReco)/eTtotalSim,1.0);
 986     if(pTtotalSim>0.0) FillHisto2D("SimPiKPPtMinusSimSmearedRecoOnly",pTtotalSim,(pTtotalSim-pTtotalReco)/pTtotalSim,1.0);
 987     if(eTtotalSim>0.0) FillHisto2D("SimPiKPEtMinusSimSmearedMultRecoOnly",nReco,(eTtotalSim-eTtotalReco)/eTtotalSim,1.0);
 988     if(pTtotalSim>0.0) FillHisto2D("SimPiKPPtMinusSimSmearedMultRecoOnly",nReco,(pTtotalSim-pTtotalReco)/pTtotalSim,1.0);
 989   }
 990   delete pID;
 991   delete strTPC;
 992   delete strITS;
 993   delete strTPCITS;
 994   return 1;
 995 }
 996 Int_t AliAnalysisHadEtMonteCarlo::AnalyseEvent(AliVEvent* ev)
 997 { // analyse MC event
 998      ResetEventValues();
 999      if(!ev){
1000             AliFatal("ERROR: Event does not exist");
1001             return 0;
1002      }
1003
1004     // Get us an mc event
1005     AliMCEvent *mcEvent = dynamic_cast<AliMCEvent*>(ev);
1006     if(!mcEvent){
1007       AliFatal("ERROR: MC Event does not exist");
1008       return 0;
1009     }
1010
1011     // Let's play with the stack!
1012     AliStack *stack = mcEvent->Stack();
1013
1014     Int_t nPrim = stack->GetNtrack();
1015
1016     Float_t fSimPiKPEtPtSmeared = 0;
1017     Float_t fSimPiKPEtEfficiencySmeared = 0;
1018     Float_t fSimPiKPEtPtCutSmearedTPC = 0;
1019     Float_t fSimPiKPEtPtCutSmearedITS = 0;
1020     Float_t fSimPiKPEtPIDSmeared = 0;
1021     Float_t fSimPiKPEtPIDSmearedNoID = 0;
1022     fSimPiKPEtShouldBeReco = 0;
1023     fSimPiKPEtShouldBeRecoPi = 0;
1024     fSimPiKPEtShouldBeRecoK = 0;
1025     fSimPiKPEtShouldBeRecoP = 0;
1026     //=================Tracks which may or may not have been reconstructed=================
1027
1028     for (Int_t iPart = 0; iPart < nPrim; iPart++)
1029     {
1030
1031       TParticle *part = stack->Particle(iPart);//This line is identified as a loss of memory by valgrind, however, the pointer still belongs to the stack, so it's the stack's problem
1032
1033         if (!part)
1034           {
1035             Printf("ERROR: Could not get particle %d", iPart);
1036             continue;
1037           }
1038         // Check if it is a primary particle
1039         if (stack->IsPhysicalPrimary(iPart)){//primaries
1040
1041           if (TMath::Abs(part->Eta()) < fHadEtReco->GetCorrections()->GetEtaCut())          {
1042
1043             Int_t pdgCode =  part->GetPDG(0)->PdgCode();
1044             bool filled = false;
1045             //Investigating smearing...
1046             //Numbers are realistic correction factors from previous studies
1047             if(fInvestigateSmearing){
1048               if(pdgCode==fgPiPlusCode ||pdgCode==fgPiMinusCode ||pdgCode==fgKPlusCode ||pdgCode==fgKMinusCode ||pdgCode==fgProtonCode ||pdgCode==fgAntiProtonCode){
1049                 //To investigate Smearing...
1050                 Float_t myet = Et(part);
1051                 fSimPiKPEt += myet;
1052                 Float_t theta = part->Theta();
1053                 Short_t charge = 1;
1054                 Float_t momentum = part->P();
1055                 //pt smearing
1056                 Float_t pSmeared = momentum *  fPtSmearer->Gaus(1,0.005);//Gaussian centered around 1
1057                 fSimPiKPEtPtSmeared += Et(pSmeared,theta,pdgCode,charge);
1058                 //Efficiency smearing
1059                 //to mock up the difference between TPC only tracks in p+p (~90% efficiency at high pT) and TPC+ITS tracks in Pb+Pb (about 70% efficiency at high pT) a factor of 7/9 was added in front
1060                 float efficiency = 7.0/9.0*2.26545*TMath::Exp(-TMath::Power(9.99977e-01/part->Pt(),7.85488e-02));//simple rough efficiency from fitting curve
1061                 if(fPtSmearer->Binomial(1,efficiency) ==1){
1062                   fSimPiKPEtEfficiencySmeared += (1.0/efficiency)*myet;
1063                 }
1064                 //pT cut smeared
1065                 if(part->Pt()>0.10){fSimPiKPEtPtCutSmearedITS +=1.00988*myet;}
1066                 if(part->Pt()>0.15){fSimPiKPEtPtCutSmearedTPC +=1.02994*myet;}
1067                 //PID smearing
1068                 fSimPiKPEtPIDSmearedNoID += 1.03018015790601458*Et(momentum,theta,fgPiPlusCode,charge);
1069                 if(part->P()<1.0){//then the particle would have been ID'd
1070                   fSimPiKPEtPIDSmeared += 1.00918051514628582*myet;
1071                 }
1072                 else{//Then it would have been assumed to be a pion
1073                   fSimPiKPEtPIDSmeared += 1.00918051514628582*Et(momentum,theta,fgPiPlusCode,charge);
1074                 }
1075               }
1076             }
1077
1078             //============Charged hadrons===================================
1079             if(pdgCode == fgPiPlusCode){
1080               float myEt = Et(part);
1081               float myEtP = Et(part,fgProtonMass);
1082               float myEtK = Et(part,fgKaonMass);
1083               if(part->Pt()>0.15) fSimPiKPEtShouldBeReco += myEt;
1084               if(part->Pt()>0.15) fSimPiKPEtShouldBeRecoPi += myEt;
1085
1086               fSimHadEt += myEt;
1087               fSimTotEt += myEt;
1088
1089
1090               if( !fRunLightweight){
1091                 FillHisto2D("EtSimulatedPiPlus",part->Pt(),part->Eta(),myEt);
1092                 FillHisto2D("EtNSimulatedPiPlus",part->Pt(),part->Eta(),1.0);
1093                 FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
1094                 FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
1095                 if(fCentBin>=0){//if a centrality bin was defined
1096                   FillHisto2D(Form("EtNSimulatedPiPlusCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1097                   FillHisto2D(Form("EtNSimulatedChargedHadronCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1098                 }
1099                 FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEt);
1100                 FillHisto2D("EtSimulatedChargedHadronAssumingProton",part->Pt(),part->Eta(),myEtP);
1101                 FillHisto2D("EtSimulatedPiPlusAssumingProton",part->Pt(),part->Eta(),myEtP);
1102                 FillHisto2D("EtSimulatedChargedHadronAssumingKaon",part->Pt(),part->Eta(),myEtK);
1103                 FillHisto2D("EtSimulatedPiPlusAssumingKaon",part->Pt(),part->Eta(),myEtK);
1104                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1105                 Short_t charge = 1;
1106                 Float_t myEtLow = Et(0.0,part->Theta(),pdgCode,charge);
1107                 Float_t myEtITS = Et(0.10/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1108                 Float_t myEtTPC = Et(0.15/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1109                 FillHisto2D("EtSimulatedChargedHadronAssumingNoPt",part->Pt(),part->Eta(),myEtLow);
1110                 FillHisto2D("EtSimulatedChargedHadronAssumingPtTPCCut",part->Pt(),part->Eta(),myEtTPC);
1111                 FillHisto2D("EtSimulatedChargedHadronAssumingPtITSCut",part->Pt(),part->Eta(),myEtITS);
1112               }
1113               filled = true;
1114             }
1115             if(pdgCode == fgPiMinusCode){
1116               float myEt = Et(part);
1117               float myEtP = Et(part,fgProtonMass);
1118               float myEtK = Et(part,fgKaonMass);
1119               if(part->Pt()>0.15) fSimPiKPEtShouldBeReco += myEt;
1120               if(part->Pt()>0.15) fSimPiKPEtShouldBeRecoPi += myEt;
1121               fSimHadEt += myEt;
1122               fSimTotEt += myEt;
1123               if( !fRunLightweight){
1124                 FillHisto2D("EtSimulatedPiMinus",part->Pt(),part->Eta(),myEt);
1125                 FillHisto2D("EtNSimulatedPiMinus",part->Pt(),part->Eta(),1.0);
1126                 FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
1127                 FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
1128                 if(fCentBin>=0){//if a centrality bin was defined
1129                   FillHisto2D(Form("EtNSimulatedPiMinusCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1130                   FillHisto2D(Form("EtNSimulatedChargedHadronCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1131                 }
1132                 FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEt);
1133                 FillHisto2D("EtSimulatedChargedHadronAssumingProton",part->Pt(),part->Eta(),myEtP);
1134                 FillHisto2D("EtSimulatedPiMinusAssumingProton",part->Pt(),part->Eta(),myEtP);
1135                 FillHisto2D("EtSimulatedChargedHadronAssumingKaon",part->Pt(),part->Eta(),myEtK);
1136                 FillHisto2D("EtSimulatedPiMinusAssumingKaon",part->Pt(),part->Eta(),myEtK);
1137                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1138                 Short_t charge = -1;
1139                 Float_t myEtLow = Et(0.0,part->Theta(),pdgCode,charge);
1140                 Float_t myEtITS = Et(0.10/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1141                 Float_t myEtTPC = Et(0.15/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1142                 FillHisto2D("EtSimulatedChargedHadronAssumingNoPt",part->Pt(),part->Eta(),myEtLow);
1143                 FillHisto2D("EtSimulatedChargedHadronAssumingPtTPCCut",part->Pt(),part->Eta(),myEtTPC);
1144                 FillHisto2D("EtSimulatedChargedHadronAssumingPtITSCut",part->Pt(),part->Eta(),myEtITS);
1145               }
1146               filled = true;
1147             }
1148             if(pdgCode == fgKPlusCode){
1149               float myEt = Et(part);
1150               float myEtPi = Et(part,fgPionMass);
1151               float myEtP = Et(part,fgProtonMass);
1152               if(part->Pt()>0.15) fSimPiKPEtShouldBeReco += myEt;
1153               if(part->Pt()>0.15) fSimPiKPEtShouldBeRecoK += myEt;
1154               fSimHadEt += myEt;
1155               fSimTotEt += myEt;
1156               if( !fRunLightweight){
1157                 FillHisto2D("EtSimulatedKPlus",part->Pt(),part->Eta(),myEt);
1158                 FillHisto2D("EtNSimulatedKPlus",part->Pt(),part->Eta(),1.0);
1159                 FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
1160                 FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
1161                 if(fCentBin>=0){//if a centrality bin was defined
1162                   FillHisto2D(Form("EtNSimulatedKPlusCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1163                   FillHisto2D(Form("EtNSimulatedChargedHadronCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1164                 }
1165                 FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
1166                 FillHisto2D("EtSimulatedKPlusAssumingPion",part->Pt(),part->Eta(),myEtPi);
1167                 FillHisto2D("EtSimulatedChargedHadronAssumingProton",part->Pt(),part->Eta(),myEtP);
1168                 FillHisto2D("EtSimulatedKPlusAssumingProton",part->Pt(),part->Eta(),myEtP);
1169                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1170                 Short_t charge = 1;
1171                 Float_t myEtLow = Et(0.0,part->Theta(),pdgCode,charge);
1172                 Float_t myEtITS = Et(0.10/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1173                 Float_t myEtTPC = Et(0.15/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1174                 FillHisto2D("EtSimulatedChargedHadronAssumingNoPt",part->Pt(),part->Eta(),myEtLow);
1175                 FillHisto2D("EtSimulatedChargedHadronAssumingPtTPCCut",part->Pt(),part->Eta(),myEtTPC);
1176                 FillHisto2D("EtSimulatedChargedHadronAssumingPtITSCut",part->Pt(),part->Eta(),myEtITS);
1177               }
1178               filled = true;
1179             }
1180             if(pdgCode == fgKMinusCode){
1181               float myEt = Et(part);
1182               float myEtPi = Et(part,fgPionMass);
1183               float myEtP = Et(part,fgProtonMass);
1184               if(part->Pt()>0.15) fSimPiKPEtShouldBeReco += myEt;
1185               if(part->Pt()>0.15) fSimPiKPEtShouldBeRecoK += myEt;
1186               fSimHadEt += myEt;
1187               fSimTotEt += myEt;
1188               if( !fRunLightweight){
1189                 FillHisto2D("EtSimulatedKMinus",part->Pt(),part->Eta(),myEt);
1190                 FillHisto2D("EtNSimulatedKMinus",part->Pt(),part->Eta(),1.0);
1191                 FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
1192                 FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
1193                 if(fCentBin>=0){//if a centrality bin was defined
1194                   FillHisto2D(Form("EtNSimulatedKMinusCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1195                   FillHisto2D(Form("EtNSimulatedChargedHadronCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1196                 }
1197                 FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
1198                 FillHisto2D("EtSimulatedKMinusAssumingPion",part->Pt(),part->Eta(),myEtPi);
1199                 FillHisto2D("EtSimulatedChargedHadronAssumingProton",part->Pt(),part->Eta(),myEtP);
1200                 FillHisto2D("EtSimulatedKMinusAssumingProton",part->Pt(),part->Eta(),myEtP);
1201                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1202                 Short_t charge = -1;
1203                 Float_t myEtLow = Et(0.0,part->Theta(),pdgCode,charge);
1204                 Float_t myEtITS = Et(0.10/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1205                 Float_t myEtTPC = Et(0.15/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1206                 FillHisto2D("EtSimulatedChargedHadronAssumingNoPt",part->Pt(),part->Eta(),myEtLow);
1207                 FillHisto2D("EtSimulatedChargedHadronAssumingPtTPCCut",part->Pt(),part->Eta(),myEtTPC);
1208                 FillHisto2D("EtSimulatedChargedHadronAssumingPtITSCut",part->Pt(),part->Eta(),myEtITS);
1209               }
1210               filled = true;
1211             }
1212             if(pdgCode == fgProtonCode){
1213               float myEt = Et(part);
1214               float myEtPi = Et(part,fgPionMass);
1215               float myEtK = Et(part,fgKaonMass);
1216               if(part->Pt()>0.15) fSimPiKPEtShouldBeReco += myEt;
1217               if(part->Pt()>0.15) fSimPiKPEtShouldBeRecoP += myEt;
1218               fSimHadEt += myEt;
1219               fSimTotEt += myEt;
1220               if( !fRunLightweight){
1221                 FillHisto2D("EtSimulatedProton",part->Pt(),part->Eta(),myEt);
1222                 FillHisto2D("EtNSimulatedProton",part->Pt(),part->Eta(),1.0);
1223                 FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
1224                 FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
1225                 if(fCentBin>=0){//if a centrality bin was defined
1226                   FillHisto2D(Form("EtNSimulatedProtonCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1227                   FillHisto2D(Form("EtNSimulatedChargedHadronCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1228                 }
1229                 FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
1230                 FillHisto2D("EtSimulatedProtonAssumingPion",part->Pt(),part->Eta(),myEtPi);
1231                 FillHisto2D("EtSimulatedChargedHadronAssumingKaon",part->Pt(),part->Eta(),myEtK);
1232                 FillHisto2D("EtSimulatedProtonAssumingKaon",part->Pt(),part->Eta(),myEtK);
1233                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1234                 Short_t charge = 1;
1235                 Float_t myEtLow = Et(0.0,part->Theta(),pdgCode,charge);
1236                 Float_t myEtITS = Et(0.10/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1237                 Float_t myEtTPC = Et(0.15/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1238                 FillHisto2D("EtSimulatedChargedHadronAssumingNoPt",part->Pt(),part->Eta(),myEtLow);
1239                 FillHisto2D("EtSimulatedChargedHadronAssumingPtTPCCut",part->Pt(),part->Eta(),myEtTPC);
1240                 FillHisto2D("EtSimulatedChargedHadronAssumingPtITSCut",part->Pt(),part->Eta(),myEtITS);
1241               }
1242               filled = true;
1243               if( !fRunLightweight){
1244                 if(fBaryonEnhancement){
1245                   float enhancement = ProtonBaryonEnhancement(part->Pt());
1246                   FillHisto2D("EtSimulatedProtonEnhanced",part->Pt(),part->Eta(),myEt*enhancement);
1247                   FillHisto2D("EtNSimulatedProtonEnhanced",part->Pt(),part->Eta(),1.0*enhancement);
1248                   FillHisto2D("EtSimulatedProtonAssumingPionEnhanced",part->Pt(),part->Eta(),myEtPi*enhancement);
1249                 }
1250               }
1251             }
1252             if(pdgCode == fgAntiProtonCode){
1253               float myEt = Et(part);
1254               float myEtPi = Et(part,fgPionMass);
1255               float myEtK = Et(part,fgKaonMass);
1256               if(part->Pt()>0.15) fSimPiKPEtShouldBeReco += myEt;
1257               if(part->Pt()>0.15) fSimPiKPEtShouldBeRecoP += myEt;
1258               fSimHadEt += myEt;
1259               fSimTotEt += myEt;
1260               if( !fRunLightweight){
1261                 FillHisto2D("EtSimulatedAntiProton",part->Pt(),part->Eta(),myEt);
1262                 FillHisto2D("EtNSimulatedAntiProton",part->Pt(),part->Eta(),1.0);
1263                 FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
1264                 FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
1265                 if(fCentBin>=0){//if a centrality bin was defined
1266                   FillHisto2D(Form("EtNSimulatedAntiProtonCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1267                   FillHisto2D(Form("EtNSimulatedChargedHadronCB%i",fCentBin),part->Pt(),part->Eta(),1.0);
1268                 }
1269                 FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
1270                 FillHisto2D("EtSimulatedAntiProtonAssumingPion",part->Pt(),part->Eta(),myEtPi);
1271                 FillHisto2D("EtSimulatedChargedHadronAssumingKaon",part->Pt(),part->Eta(),myEtK);
1272                 FillHisto2D("EtSimulatedAntiProtonAssumingKaon",part->Pt(),part->Eta(),myEtK);
1273                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1274                 Short_t charge = -1;
1275                 Float_t myEtLow = Et(0.0,part->Theta(),pdgCode,charge);
1276                 Float_t myEtITS = Et(0.10/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1277                 Float_t myEtTPC = Et(0.15/TMath::Sin(part->Theta()),part->Theta(),pdgCode,charge);
1278                 FillHisto2D("EtSimulatedChargedHadronAssumingNoPt",part->Pt(),part->Eta(),myEtLow);
1279                 FillHisto2D("EtSimulatedChargedHadronAssumingPtTPCCut",part->Pt(),part->Eta(),myEtTPC);
1280                 FillHisto2D("EtSimulatedChargedHadronAssumingPtITSCut",part->Pt(),part->Eta(),myEtITS);
1281               }
1282               filled = true;
1283               if( !fRunLightweight){
1284                 if(fBaryonEnhancement){
1285                   float enhancement = ProtonBaryonEnhancement(part->Pt());
1286                   FillHisto2D("EtSimulatedAntiProtonEnhanced",part->Pt(),part->Eta(),myEt*enhancement);
1287                   FillHisto2D("EtNSimulatedAntiProtonEnhanced",part->Pt(),part->Eta(),1.0*enhancement);
1288                   FillHisto2D("EtSimulatedAntiProtonAssumingPionEnhanced",part->Pt(),part->Eta(),myEtPi*enhancement);
1289                 }
1290               }
1291             }
1292             //============Other hadrons===================================
1293
1294             if(pdgCode == fgNeutronCode){
1295               float myEt = Et(part);
1296               fSimHadEt += myEt;
1297               fSimTotEt += myEt;
1298               if( !fRunLightweight){
1299                 FillHisto2D("EtSimulatedNeutron",part->Pt(),part->Eta(),myEt);
1300                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1301               }
1302               filled = true;
1303             }
1304             if(pdgCode == fgAntiNeutronCode){
1305               float myEt = Et(part);
1306               fSimHadEt += myEt;
1307               fSimTotEt += myEt;
1308               if( !fRunLightweight){
1309                 FillHisto2D("EtSimulatedAntiNeutron",part->Pt(),part->Eta(),myEt);
1310                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1311               }
1312               filled = true;
1313             }
1314             if(pdgCode == fgLambdaCode){
1315               float myEt = Et(part);
1316               fSimHadEt += myEt;
1317               fSimTotEt += myEt;
1318               if( !fRunLightweight){
1319                 FillHisto2D("EtSimulatedLambda",part->Pt(),part->Eta(),myEt);
1320                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1321                 Float_t weight = LambdaWeight(part->Pt());
1322                 if(fBaryonEnhancement){
1323                   float enhancement = ProtonBaryonEnhancement(part->Pt());
1324                   weight = weight*enhancement;
1325                 }
1326                 FillHisto2D("EtSimulatedLambdaReweighted",part->Pt(),part->Eta(),myEt*weight);
1327                 Int_t ndaughters = part->GetNDaughters();
1328                 for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
1329                   Int_t daughterindex = part->GetDaughter(idaughter);
1330                   if(daughterindex<0 || daughterindex>1e5) continue;
1331                   TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
1332                   if(daughter){
1333                     if(daughter->GetPDG(0)){
1334                       Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
1335                       if(daughtercode==fgPiMinusCode || daughtercode==fgProtonCode){
1336                         myEt = Et(daughter);
1337                         FillHisto2D("EtSimulatedLambdaDaughters",daughter->Pt(),daughter->Eta(),myEt);
1338                         FillHisto2D("EtSimulatedLambdaDaughtersReweighted",daughter->Pt(),daughter->Eta(),myEt*weight);
1339                       }
1340                     }
1341                     else{
1342                     }
1343                   }
1344                 }
1345               }
1346               filled = true;
1347             }
1348             if(pdgCode == fgAntiLambdaCode){
1349               float myEt = Et(part);
1350               fSimHadEt += myEt;
1351               fSimTotEt += myEt;
1352               if( !fRunLightweight){
1353                 FillHisto2D("EtSimulatedAntiLambda",part->Pt(),part->Eta(),myEt);
1354                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1355                 Float_t weight = AntiLambdaWeight(part->Pt());
1356                 if(fBaryonEnhancement){
1357                   float enhancement = ProtonBaryonEnhancement(part->Pt());
1358                   weight = weight*enhancement;
1359                 }
1360                 FillHisto2D("EtSimulatedAntiLambdaReweighted",part->Pt(),part->Eta(),myEt*weight);
1361                 Int_t ndaughters = part->GetNDaughters();
1362                 for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
1363                   Int_t daughterindex = part->GetDaughter(idaughter);
1364                   if(daughterindex<0 || daughterindex>1e5) continue;
1365                   TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
1366                   if(daughter){
1367                     if(daughter->GetPDG(0)){
1368                       Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
1369                       if(daughtercode==fgPiPlusCode || daughtercode==fgAntiProtonCode){
1370                         myEt = Et(daughter);
1371                         FillHisto2D("EtSimulatedAntiLambdaDaughters",daughter->Pt(),daughter->Eta(),myEt);
1372                         FillHisto2D("EtSimulatedAntiLambdaDaughtersReweighted",daughter->Pt(),daughter->Eta(),myEt*weight);
1373                       }
1374                     }
1375                   }
1376                 }
1377               }
1378               filled = true;
1379             }
1380             if(pdgCode == fgK0SCode){
1381               float myEt = Et(part);
1382               fSimHadEt += myEt;
1383               fSimTotEt += myEt;
1384               if( !fRunLightweight){
1385                 FillHisto2D("EtSimulatedK0S",part->Pt(),part->Eta(),myEt);
1386                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1387                 Float_t weight = K0Weight(part->Pt());
1388                 FillHisto2D("EtSimulatedK0SReweighted",part->Pt(),part->Eta(),myEt*weight);
1389                 Int_t ndaughters = part->GetNDaughters();
1390                 for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
1391                   Int_t daughterindex = part->GetDaughter(idaughter);
1392                   if(daughterindex<0 || daughterindex>1e5) continue;
1393                   TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
1394                   if(daughter){
1395                     if(daughter->GetPDG(0)){
1396
1397                       Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
1398                       if(daughtercode==fgPiMinusCode || daughtercode==fgPiPlusCode){
1399                         myEt = Et(daughter);
1400                         FillHisto2D("EtSimulatedK0SDaughters",daughter->Pt(),daughter->Eta(),myEt);
1401                         FillHisto2D("EtSimulatedK0SDaughtersReweighted",daughter->Pt(),daughter->Eta(),myEt*weight);
1402                       }
1403                     }
1404                   }
1405                 }
1406               }
1407               filled = true;
1408             }
1409             if(pdgCode == fgK0LCode){
1410               float myEt = Et(part);
1411               fSimHadEt += myEt;
1412               fSimTotEt += myEt;
1413               if( !fRunLightweight){
1414                 FillHisto2D("EtSimulatedK0L",part->Pt(),part->Eta(),myEt);
1415                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1416                 Float_t weight = K0Weight(part->Pt());
1417                 FillHisto2D("EtSimulatedK0LReweighted",part->Pt(),part->Eta(),myEt*weight);
1418               }
1419               filled = true;
1420             }
1421             if(pdgCode == fgOmegaCode){
1422               float myEt = Et(part);
1423               fSimHadEt += myEt;
1424               fSimTotEt += myEt;
1425               if( !fRunLightweight){
1426                 FillHisto2D("EtSimulatedOmega",part->Pt(),part->Eta(),myEt);
1427                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1428                 Int_t ndaughters = part->GetNDaughters();
1429                 for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
1430                   Int_t daughterindex = part->GetDaughter(idaughter);
1431                   if(daughterindex<0 || daughterindex>1e5) continue;
1432                   TParticle *daughter = stack->Particle(daughterindex);
1433                   if(daughter){
1434                     if(daughter->GetPDG(0)){
1435
1436                       Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
1437                       if(daughtercode==fgPiPlusCode || daughtercode==fgProtonCode || daughtercode==fgKMinusCode){
1438                         myEt = Et(daughter);
1439                         FillHisto2D("EtSimulatedOmegaDaughters",daughter->Pt(),daughter->Eta(),myEt);
1440                       }
1441                     }
1442                   }
1443                 }
1444               }
1445               filled = true;
1446             }
1447             if(pdgCode == fgAntiOmegaCode){
1448               float myEt = Et(part);
1449               fSimHadEt += myEt;
1450               fSimTotEt += myEt;
1451               if( !fRunLightweight){
1452                 FillHisto2D("EtSimulatedOmega",part->Pt(),part->Eta(),myEt);
1453                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1454                 Int_t ndaughters = part->GetNDaughters();
1455                 for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
1456                   Int_t daughterindex = part->GetDaughter(idaughter);
1457                   if(daughterindex<0 || daughterindex>1e5) continue;
1458                   TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
1459                   if(daughter){
1460                     if(daughter->GetPDG(0)){
1461                       Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
1462                       if(daughtercode==fgPiMinusCode || daughtercode==fgAntiProtonCode || daughtercode==fgKPlusCode){
1463                         myEt = Et(daughter);
1464                         FillHisto2D("EtSimulatedAntiOmegaDaughters",daughter->Pt(),daughter->Eta(),myEt);
1465                       }
1466                     }
1467                   }
1468                 }
1469               }
1470               filled = true;
1471             }
1472             //There are two codes for Sigmas
1473             if(pdgCode == fgSigmaCode || pdgCode == -3222){
1474               float myEt = Et(part);
1475               fSimHadEt += myEt;
1476               fSimTotEt += myEt;
1477               if( !fRunLightweight){
1478                 FillHisto2D("EtSimulatedSigma",part->Pt(),part->Eta(),myEt);
1479                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1480               }
1481               filled = true;
1482             }
1483             if(pdgCode == fgAntiSigmaCode || pdgCode == 3222){
1484               float myEt = Et(part);
1485               fSimHadEt += myEt;
1486               fSimTotEt += myEt;
1487               if( !fRunLightweight){
1488                 FillHisto2D("EtSimulatedAntiSigma",part->Pt(),part->Eta(),myEt);
1489                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1490               }
1491               filled = true;
1492             }
1493             if(pdgCode == fgXiCode){
1494               float myEt = Et(part);
1495               fSimHadEt += myEt;
1496               fSimTotEt += myEt;
1497               if( !fRunLightweight){
1498                 FillHisto2D("EtSimulatedXi",part->Pt(),part->Eta(),myEt);
1499                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1500                 Int_t ndaughters = part->GetNDaughters();
1501                 for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
1502                   Int_t daughterindex = part->GetDaughter(idaughter);
1503                   if(daughterindex<0 || daughterindex>1e5 || daughterindex>1e5) continue;
1504                   TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
1505                   if(daughter){
1506                     if(daughter->GetPDG(0)){
1507
1508                       Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
1509                       if(daughtercode==fgPiPlusCode || daughtercode==fgProtonCode || daughtercode==fgPiMinusCode){
1510                         myEt = Et(daughter);
1511                         FillHisto2D("EtSimulatedXiDaughters",daughter->Pt(),daughter->Eta(),myEt);
1512                       }
1513                     }
1514                   }
1515                 }
1516               }
1517               filled = true;
1518             }
1519             if(pdgCode == fgAntiXiCode){
1520               float myEt = Et(part);
1521               fSimHadEt += myEt;
1522               fSimTotEt += myEt;
1523               if( !fRunLightweight){
1524                 FillHisto2D("EtSimulatedAntiXi",part->Pt(),part->Eta(),myEt);
1525                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1526                 Int_t ndaughters = part->GetNDaughters();
1527                 for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
1528                   Int_t daughterindex = part->GetDaughter(idaughter);
1529                   if(daughterindex<0 || daughterindex>1e5) continue;
1530                   TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
1531                   if(daughter){
1532                     if(daughter->GetPDG(0)){
1533                       Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
1534                       if(daughtercode==fgPiPlusCode || daughtercode==fgAntiProtonCode || daughtercode==fgPiMinusCode){
1535                         myEt = Et(daughter);
1536                         FillHisto2D("EtSimulatedAntiXiDaughters",daughter->Pt(),daughter->Eta(),myEt);
1537                       }
1538                     }
1539                   }
1540                 }
1541               }
1542               filled = true;
1543             }
1544             if(pdgCode == fgXi0Code){
1545               float myEt = Et(part);
1546               fSimHadEt += myEt;
1547               fSimTotEt += myEt;
1548               if( !fRunLightweight){
1549                 FillHisto2D("EtSimulatedXi0",part->Pt(),part->Eta(),myEt);
1550                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1551               }
1552               filled = true;
1553             }
1554             if(pdgCode == fgAntiXi0Code){
1555               float myEt = Et(part);
1556               fSimHadEt += myEt;
1557               fSimTotEt += myEt;
1558               if( !fRunLightweight){
1559                 FillHisto2D("EtSimulatedAntiXi0",part->Pt(),part->Eta(),myEt);
1560                 FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
1561               }
1562               filled = true;
1563             }
1564             //============electrons===================================
1565
1566             if(pdgCode == fgEPlusCode){
1567               float myEt = Et(part);
1568               fSimTotEt += myEt;
1569               if( !fRunLightweight){
1570                 FillHisto2D("EtSimulatedEPlus",part->Pt(),part->Eta(),myEt);
1571               }
1572               filled = true;
1573             }
1574             if(pdgCode == fgEMinusCode){
1575               float myEt = Et(part);
1576               fSimTotEt += myEt;
1577               if( !fRunLightweight){
1578                 FillHisto2D("EtSimulatedEMinus",part->Pt(),part->Eta(),myEt);
1579               }
1580               filled = true;
1581             }
1582             //============neutrals===================================
1583             if(pdgCode == fgGammaCode){
1584               TParticle *mom = NULL;
1585               Int_t pdgCodeMom = -99999999;
1586               float momEta = -30;
1587               float mompT = -5;
1588               if(part->GetFirstMother()){
1589                 mom = stack->Particle(part->GetFirstMother());
1590                 pdgCodeMom =  mom->GetPDG(0)->PdgCode();
1591                 momEta = mom->Eta();
1592                 mompT = mom->Pt();
1593               }
1594               //We want to separate the gammas by pi0, eta, omega0 but we don't want to double count energy so we get the et from the gamma daughter
1595               if(pdgCodeMom == fgEtaCode){
1596                 float myEt = Et(part);
1597                 fSimTotEt += myEt;
1598                 if( !fRunLightweight){
1599                   FillHisto2D("EtSimulatedEta",mompT,momEta,myEt);
1600                 }
1601                 filled = true;
1602               }
1603               if(pdgCodeMom == fgPi0Code){
1604                 float myEt = Et(part);
1605                 fSimTotEt += myEt;
1606                 if( !fRunLightweight){
1607                   FillHisto2D("EtSimulatedPi0",mompT,momEta,myEt);
1608                 }
1609                 filled = true;
1610               }
1611               if(pdgCodeMom == fgOmega0Code){
1612                 float myEt = Et(part);
1613                 fSimTotEt += myEt;
1614                 if( !fRunLightweight){
1615                   FillHisto2D("EtSimulatedOmega0",mompT,momEta,myEt);
1616                 }
1617                 filled = true;
1618               }
1619               if(!filled){
1620                 float myEt = Et(part);
1621                 fSimTotEt += myEt;
1622                 if( !fRunLightweight){
1623                   FillHisto2D("EtSimulatedGamma",part->Pt(),part->Eta(),myEt);
1624                 }
1625                 filled = true;
1626               }
1627             }
1628             if(pdgCode == fgEtaCode){
1629               float myEt = Et(part);
1630               fSimTotEt += myEt;
1631               if( !fRunLightweight){
1632                 FillHisto2D("EtSimulatedEta",part->Pt(),part->Eta(),myEt);
1633               }
1634               filled = true;
1635             }
1636             if(pdgCode == fgPi0Code){
1637               float myEt = Et(part);
1638               fSimTotEt += myEt;
1639               if( !fRunLightweight){
1640                 FillHisto2D("EtSimulatedPi0",part->Pt(),part->Eta(),myEt);
1641               }
1642               filled = true;
1643             }
1644             if(pdgCode == fgOmega0Code){
1645               float myEt = Et(part);
1646               fSimTotEt += myEt;
1647               if( !fRunLightweight){
1648                 FillHisto2D("EtSimulatedOmega0",part->Pt(),part->Eta(),myEt);
1649               }
1650               filled = true;
1651             }
1652           }
1653         }
1654     }
1655
1656     FillHisto1D("SimTotEt",fSimTotEt,1.0);
1657     FillHisto1D("SimHadEt",fSimHadEt,1.0);
1658     FillHisto1D("SimPiKPEt",fSimPiKPEt,1.0);
1659     if(AliPWG0Helper::GetEventProcessType(mcEvent->Header()) == AliPWG0Helper::kND){
1660       FillHisto1D("SimHadEtND",fSimHadEt,1.0);
1661       FillHisto1D("SimTotEtND",fSimHadEt,1.0);
1662       FillHisto1D("NEventsND",0.5,1);
1663       FillHisto1D("SimPiKPEtND",fSimPiKPEt,1.0);
1664     }
1665     if(AliPWG0Helper::GetEventProcessType(mcEvent->Header()) == AliPWG0Helper::kSD){
1666       FillHisto1D("SimHadEtSD",fSimHadEt,1.0);
1667       FillHisto1D("SimTotEtSD",fSimHadEt,1.0);
1668       FillHisto1D("NEventsSD",0.5,1);
1669       FillHisto1D("SimPiKPEtSD",fSimPiKPEt,1.0);
1670     }
1671     if(AliPWG0Helper::GetEventProcessType(mcEvent->Header()) == AliPWG0Helper::kDD){
1672       FillHisto1D("SimHadEtDD",fSimHadEt,1.0);
1673       FillHisto1D("SimTotEtDD",fSimHadEt,1.0);
1674       FillHisto1D("NEventsDD",0.5,1);
1675       FillHisto1D("SimPiKPEtSD",fSimPiKPEt,1.0);
1676     }
1677     if(fCentBin != -1){//if we have Pb+Pb and a centrality bin was found
1678       if(fSimTotEt>0.0) FillHisto1D(Form("SimTotEtCB%i",fCentBin),fSimTotEt,1.0);
1679       if(fSimHadEt>0.0) FillHisto1D(Form("SimHadEtCB%i",fCentBin),fSimHadEt,1.0);
1680       if(fSimPiKPEt>0.0)FillHisto1D(Form("SimPiKPEtCB%i",fCentBin),fSimPiKPEt,1.0);
1681     }
1682
1683     if(fInvestigateSmearing && !fRunLightweight){
1684       //Smearing histograms
1685       if(fSimPiKPEt>0.0) FillHisto2D("SimPiKPEtMinusSimPtSmeared",fSimPiKPEt,(fSimPiKPEt-fSimPiKPEtPtSmeared)/fSimPiKPEt,1.0);
1686       FillHisto1D("SimPiKPEtPtSmeared",fSimPiKPEtPtSmeared,1.0);
1687       if(fSimPiKPEt>0.0) FillHisto2D("SimPiKPEtMinusSimEfficiencySmeared",fSimPiKPEt,(fSimPiKPEt-fSimPiKPEtEfficiencySmeared)/fSimPiKPEt,1.0);
1688       FillHisto1D("SimPiKPEtEfficiencySmeared",fSimPiKPEtEfficiencySmeared,1.0);
1689       if(fSimPiKPEt>0.0) FillHisto2D("SimPiKPEtMinusSimPtCutSmearedTPC",fSimPiKPEt,(fSimPiKPEt-fSimPiKPEtPtCutSmearedTPC)/fSimPiKPEt,1.0);
1690       FillHisto1D("SimPiKPEtPtCutSmearedTPC",fSimPiKPEtPtCutSmearedTPC,1.0);
1691       if(fSimPiKPEt>0.0) FillHisto2D("SimPiKPEtMinusSimPtCutSmearedITS",fSimPiKPEt,(fSimPiKPEt-fSimPiKPEtPtCutSmearedITS)/fSimPiKPEt,1.0);
1692       FillHisto1D("SimPiKPEtPtCutSmearedITS",fSimPiKPEtPtCutSmearedTPC,1.0);
1693       if(fSimPiKPEt>0.0) FillHisto2D("SimPiKPEtMinusSimPIDSmeared",fSimPiKPEt,(fSimPiKPEt-fSimPiKPEtPIDSmeared)/fSimPiKPEt,1.0);
1694       FillHisto1D("SimPiKPEtPIDSmeared",fSimPiKPEtPIDSmeared,1.0);
1695       if(fSimPiKPEt>0.0) FillHisto2D("SimPiKPEtMinusSimPIDSmearedNoID",fSimPiKPEt,(fSimPiKPEt-fSimPiKPEtPIDSmearedNoID)/fSimPiKPEt,1.0);
1696       FillHisto1D("SimPiKPEtPIDSmearedNoID",fSimPiKPEtPIDSmearedNoID,1.0);
1697     }
1698     return 1;
1699
1700 }
1701
1702 void AliAnalysisHadEtMonteCarlo::Init()
1703 { // Init
1704     AliAnalysisHadEt::Init();
1705     if(!fPtSmearer) fPtSmearer = new TRandom();
1706 }
1707 void AliAnalysisHadEtMonteCarlo::CreateHistograms(){
1708   //for simulated Et only (no reconstruction)
1709   if( !fRunLightweight){
1710     CreateEtaPtHisto2D(TString("EtSimulatedPiPlus"),TString("Simulated E_{T} from #pi^{+}"));
1711     CreateEtaPtHisto2D("EtSimulatedPiMinus","Simulated E_{T} from #pi^{-}");
1712     CreateEtaPtHisto2D("EtSimulatedKPlus","Simulated E_{T} from K^{+}");
1713     CreateEtaPtHisto2D("EtSimulatedKMinus","Simulated E_{T} from K^{-}");
1714     CreateEtaPtHisto2D("EtSimulatedProton","Simulated E_{T} from p");
1715     CreateEtaPtHisto2D("EtSimulatedAntiProton","Simulated E_{T} from #bar{p}");//Both baryon enhancement and strangeness rescaling
1716     if(fBaryonEnhancement){
1717       CreateEtaPtHisto2D("EtSimulatedProtonEnhanced","Simulated E_{T} from p");
1718       CreateEtaPtHisto2D("EtSimulatedAntiProtonEnhanced","Simulated E_{T} from #bar{p}");
1719     }
1720     CreateEtaPtHisto2D("EtSimulatedChargedHadron","Simulated E_{T} from charged hadrons");
1721     CreateEtaPtHisto2D("EtNSimulatedPiPlus","Number of Simulated #pi^{+}");
1722     CreateEtaPtHisto2D("EtNSimulatedPiMinus","Number of simulated #pi^{-}");
1723     CreateEtaPtHisto2D("EtNSimulatedKPlus","Number of simulated K^{+}");
1724     CreateEtaPtHisto2D("EtNSimulatedKMinus","Number of simulated K^{-}");
1725     CreateEtaPtHisto2D("EtNSimulatedProton","Number of simulated p");
1726     CreateEtaPtHisto2D("EtNSimulatedAntiProton","Number of simulated #bar{p}");
1727     if(fBaryonEnhancement){
1728       CreateEtaPtHisto2D("EtNSimulatedProtonEnhanced","Number of simulated p");
1729       CreateEtaPtHisto2D("EtNSimulatedAntiProtonEnhanced","Number of simulated #bar{p}");
1730     }
1731     CreateEtaPtHisto2D("EtNSimulatedChargedHadron","Number of simulated charged hadrons");
1732     if(fDataSet==20100){//If this is Pb+Pb
1733       Int_t width = 5;
1734       if(fNCentBins<21) width = 10;
1735       for(Int_t i=0;i<fNCentBins;i++){
1736         CreateEtaPtHisto2D(Form("EtNSimulatedPiPlusCB%i",i),Form("Number of Simulated #pi^{+} for %i-%i central",i*width,(i+1)*width));
1737         CreateEtaPtHisto2D(Form("EtNSimulatedPiMinusCB%i",i),Form("Number of simulated #pi^{-} for %i-%i central",i*width,(i+1)*width));
1738         CreateEtaPtHisto2D(Form("EtNSimulatedKPlusCB%i",i),Form("Number of simulated K^{+} for %i-%i central",i*width,(i+1)*width));
1739         CreateEtaPtHisto2D(Form("EtNSimulatedKMinusCB%i",i),Form("Number of simulated K^{-} for %i-%i central",i*width,(i+1)*width));
1740         CreateEtaPtHisto2D(Form("EtNSimulatedProtonCB%i",i),Form("Number of simulated p for %i-%i central",i*width,(i+1)*width));
1741         CreateEtaPtHisto2D(Form("EtNSimulatedAntiProtonCB%i",i),Form("Number of simulated #bar{p} for %i-%i central",i*width,(i+1)*width));
1742         CreateEtaPtHisto2D(Form("EtNSimulatedChargedHadronCB%i",i),Form("Number of simulated charged hadrons for %i-%i central",i*width,(i+1)*width));
1743       }
1744     }
1745     CreateEtaPtHisto2D("EtSimulatedChargedHadronAssumingNoPt","Simulated E_{T} from charged hadrons assuming p_{T}=0");
1746     CreateEtaPtHisto2D("EtSimulatedChargedHadronAssumingPtTPCCut","Simulated E_{T} from charged hadrons assuming p_{T}=0.15");
1747     CreateEtaPtHisto2D("EtSimulatedChargedHadronAssumingPtITSCut","Simulated E_{T} from charged hadrons assuming p_{T}=0.10");
1748
1749     CreateEtaPtHisto2D("EtSimulatedChargedHadronAssumingPion","Simulated E_{T} from charged hadrons assuming they are all pions");
1750     CreateEtaPtHisto2D("EtSimulatedChargedHadronAssumingProton","Simulated E_{T} from charged hadrons assuming they are all pions");
1751     CreateEtaPtHisto2D("EtSimulatedChargedHadronAssumingKaon","Simulated E_{T} from charged hadrons assuming they are all pions");
1752     CreateEtaPtHisto2D("EtSimulatedKPlusAssumingPion","Simulated E_{T} from K^{+} assuming #pi mass");
1753     CreateEtaPtHisto2D("EtSimulatedKMinusAssumingPion","Simulated E_{T} from K^{-} assuming #pi mass");
1754     CreateEtaPtHisto2D("EtSimulatedProtonAssumingPion","Simulated E_{T} from p assuming #pi mass");
1755     CreateEtaPtHisto2D("EtSimulatedAntiProtonAssumingPion","Simulated E_{T} from #bar{p} assuming #pi mass");
1756     CreateEtaPtHisto2D("EtSimulatedKPlusAssumingProton","Simulated E_{T} from K^{+} assuming #pi mass");
1757     CreateEtaPtHisto2D("EtSimulatedKMinusAssumingProton","Simulated E_{T} from K^{-} assuming #pi mass");
1758     CreateEtaPtHisto2D("EtSimulatedPiPlusAssumingProton","Simulated E_{T} from p assuming #pi mass");
1759     CreateEtaPtHisto2D("EtSimulatedPiMinusAssumingProton","Simulated E_{T} from #bar{p} assuming #pi mass");
1760     CreateEtaPtHisto2D("EtSimulatedPiPlusAssumingKaon","Simulated E_{T} from K^{+} assuming #pi mass");
1761     CreateEtaPtHisto2D("EtSimulatedPiMinusAssumingKaon","Simulated E_{T} from K^{-} assuming #pi mass");
1762     CreateEtaPtHisto2D("EtSimulatedProtonAssumingKaon","Simulated E_{T} from p assuming #pi mass");
1763     CreateEtaPtHisto2D("EtSimulatedAntiProtonAssumingKaon","Simulated E_{T} from #bar{p} assuming #pi mass");
1764     if(fBaryonEnhancement){
1765       CreateEtaPtHisto2D("EtSimulatedProtonAssumingPionEnhanced","Simulated E_{T} from p assuming #pi mass");
1766       CreateEtaPtHisto2D("EtSimulatedAntiProtonAssumingPionEnhanced","Simulated E_{T} from #bar{p} assuming #pi mass");
1767     }
1768
1769     CreateEtaPtHisto2D("EtSimulatedLambda","Simulated E_{T} from #Lambda");
1770     CreateEtaPtHisto2D("EtSimulatedAntiLambda","Simulated E_{T} from #bar{#Lambda}");
1771     CreateEtaPtHisto2D("EtSimulatedK0S","Simulated E_{T} from K^{0}_{S}");
1772     CreateEtaPtHisto2D("EtSimulatedK0L","Simulated E_{T} from K^{0}_{L}");
1773     CreateEtaPtHisto2D("EtSimulatedLambdaReweighted","Simulated E_{T} from #Lambda");//These will also be used for baryon enhancement
1774     CreateEtaPtHisto2D("EtSimulatedAntiLambdaReweighted","Simulated E_{T} from #bar{#Lambda}");
1775     CreateEtaPtHisto2D("EtSimulatedK0SReweighted","Simulated E_{T} from K^{0}_{S}");
1776     CreateEtaPtHisto2D("EtSimulatedK0LReweighted","Simulated E_{T} from K^{0}_{L}");
1777     CreateEtaPtHisto2D("EtSimulatedNeutron","Simulated E_{T} from neutrons");
1778     CreateEtaPtHisto2D("EtSimulatedAntiNeutron","Simulated E_{T} from #bar{n}");
1779     CreateEtaPtHisto2D("EtSimulatedEPlus","Simulated E_{T} from e^{+}");
1780     CreateEtaPtHisto2D("EtSimulatedEMinus","Simulated E_{T} from e^{-}");
1781     CreateEtaPtHisto2D("EtSimulatedOmega","Simulated E_{T} from #Omega^{-}");
1782     CreateEtaPtHisto2D("EtSimulatedAntiOmega","Simulated E_{T} from #Omega^{+}");
1783     CreateEtaPtHisto2D("EtSimulatedXi","Simulated E_{T} from #Xi^{-}");
1784     CreateEtaPtHisto2D("EtSimulatedAntiXi","Simulated E_{T} from #Xi^{+}");
1785     CreateEtaPtHisto2D("EtSimulatedSigma","Simulated E_{T} from #Xi^{-}");
1786     CreateEtaPtHisto2D("EtSimulatedAntiSigma","Simulated E_{T} from #Xi^{+}");
1787     CreateEtaPtHisto2D("EtSimulatedXi0","Simulated E_{T} from #Xi^{0}");
1788     CreateEtaPtHisto2D("EtSimulatedAntiXi0","Simulated E_{T} from #Xi^{0}");
1789     CreateEtaPtHisto2D("EtSimulatedAllHadron","Simulated E_{T} from all hadrons");
1790
1791
1792     CreateEtaPtHisto2D("EtSimulatedLambdaDaughters","Simulated E_{T} from #Lambda Daughters");
1793     CreateEtaPtHisto2D("EtSimulatedAntiLambdaDaughters","Simulated E_{T} from #bar{#Lambda} Daughters");
1794     CreateEtaPtHisto2D("EtSimulatedK0SDaughters","Simulated E_{T} from K^{0}_{S} Daughters");
1795     CreateEtaPtHisto2D("EtSimulatedLambdaDaughtersReweighted","Simulated E_{T} from #Lambda Daughters");
1796     CreateEtaPtHisto2D("EtSimulatedAntiLambdaDaughtersReweighted","Simulated E_{T} from #bar{#Lambda} Daughters");
1797     CreateEtaPtHisto2D("EtSimulatedK0SDaughtersReweighted","Simulated E_{T} from K^{0}_{S} Daughters");
1798     CreateEtaPtHisto2D("EtSimulatedOmegaDaughters","Simulated E_{T} from #Omega^{-} Daughters");
1799     CreateEtaPtHisto2D("EtSimulatedAntiOmegaDaughters","Simulated E_{T} from #Omega^{+} Daughters");
1800     CreateEtaPtHisto2D("EtSimulatedXiDaughters","Simulated E_{T} from #Xi^{-} Daughters");
1801     CreateEtaPtHisto2D("EtSimulatedAntiXiDaughters","Simulated E_{T} from #Xi^{+} Daughters");
1802
1803
1804     CreateEtaPtHisto2D("EtSimulatedGamma","Simulated E_{T} from #gamma");
1805     CreateEtaPtHisto2D("EtSimulatedEta","Simulated E_{T} from #eta");
1806     CreateEtaPtHisto2D("EtSimulatedPi0","Simulated E_{T} from #pi^{0}");
1807     CreateEtaPtHisto2D("EtSimulatedOmega0","Simulated E_{T} from #omega");
1808   }
1809   TString *strTPC = new TString("TPC");
1810   TString *strITS = new TString("ITS");
1811   TString *strTPCITS = new TString("TPCITS");
1812   Int_t lastcutset = 1;
1813   if(fRequireITSHits) lastcutset = 2;
1814   for(Int_t i=0;i<=lastcutset;i++){
1815     TString *cutName = NULL;
1816     Float_t maxPtdEdx = 10;
1817     Float_t mindEdx = 35;
1818     Float_t maxdEdx = 150.0;
1819     switch(i){
1820     case 0:
1821       cutName = strTPC;
1822       break;
1823     case 1:
1824       cutName = strITS;
1825       maxPtdEdx = 5;
1826       maxdEdx = 500.0;
1827       break;
1828     case 2:
1829       cutName = strTPCITS;
1830       break;
1831     default:
1832       cerr<<"Error:  cannot make histograms!"<<endl;
1833       return;
1834     }
1835
1836     if( !fRunLightweight){
1837       CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedPiPlus",cutName->Data()),"Reconstructed E_{T} from identified #pi^{+}");
1838       CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedPiMinus",cutName->Data()),"Reconstructed E_{T} from identified #pi^{-}");
1839       CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedKPlus",cutName->Data()),"Reconstructed E_{T} from identified K^{+}");
1840       CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedEMinus",cutName->Data()),"Reconstructed E_{T} from identified e^{-}");
1841       CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedEPlus",cutName->Data()),"Reconstructed E_{T} from identified e^{+}");
1842       CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedKMinus",cutName->Data()),"Reconstructed E_{T} from identified K^{-}");
1843       CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedProton",cutName->Data()),"Reconstructed E_{T} from identified p");
1844       CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedAntiProton",cutName->Data()),"Reconstructed E_{T} from identified #bar{p}");
1845       if(fBaryonEnhancement){
1846         CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedProtonEnhanced",cutName->Data()),"Reconstructed E_{T} from identified p");
1847         CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedAntiProtonEnhanced",cutName->Data()),"Reconstructed E_{T} from identified #bar{p}");
1848       }
1849       CreateEtaPtHisto2D(Form("EtNReconstructed%sUnidentified",cutName->Data()),"Number of Reconstructed unidentified particles");
1850       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedAssumingPion",cutName->Data()),"Reconstructed E_{T} from unidentified particles assuming pion mass");
1851       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedAssumingProton",cutName->Data()),"Reconstructed E_{T} from unidentified particles assuming pion mass");
1852       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedAssumingKaon",cutName->Data()),"Reconstructed E_{T} from unidentified particles assuming pion mass");
1853
1854       CreateEtaPtHisto2D(Form("EtNReconstructed%sUnidentifiedKaon",cutName->Data()),"Number of Reconstructed unidentified kaons particles");
1855       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedKaonAssumingPion",cutName->Data()),"Reconstructed E_{T} from unidentified kaons particles assuming pion mass");
1856       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedKaonAssumingProton",cutName->Data()),"Reconstructed E_{T} from unidentified kaons particles assuming pion mass");
1857       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedKaonAssumingKaon",cutName->Data()),"Reconstructed E_{T} from unidentified kaons particles assuming pion mass");
1858       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedKaon",cutName->Data()),"Reconstructed E_{T} from unidentified kaons particles assuming kaon mass");
1859       CreateEtaPtHisto2D(Form("EtNReconstructed%sUnidentifiedProton",cutName->Data()),"Number of Reconstructed unidentified proton particles");
1860       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedProtonAssumingPion",cutName->Data()),"Reconstructed E_{T} from unidentified proton particles assuming pion mass");
1861       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedProtonAssumingKaon",cutName->Data()),"Reconstructed E_{T} from unidentified proton particles assuming pion mass");
1862       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedPionAssumingKaon",cutName->Data()),"Reconstructed E_{T} from unidentified kaons particles assuming pion mass");
1863       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedProtonAssumingKaon",cutName->Data()),"Reconstructed E_{T} from unidentified proton particles assuming pion mass");
1864       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedPionAssumingProton",cutName->Data()),"Reconstructed E_{T} from unidentified kaons particles assuming pion mass");
1865       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedProtonAssumingProton",cutName->Data()),"Reconstructed E_{T} from unidentified proton particles assuming pion mass");
1866       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedProton",cutName->Data()),"Reconstructed E_{T} from unidentified proton particles assuming proton mass");
1867       if(fBaryonEnhancement){
1868         CreateEtaPtHisto2D(Form("EtNReconstructed%sUnidentifiedProtonEnhanced",cutName->Data()),"Number of Reconstructed unidentified proton particles");
1869         CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedProtonAssumingPionEnhanced",cutName->Data()),"Reconstructed E_{T} from unidentified proton particles assuming pion mass");
1870         CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedProtonEnhanced",cutName->Data()),"Reconstructed E_{T} from unidentified proton particles assuming proton mass");
1871       }
1872       CreateEtaPtHisto2D(Form("EtNReconstructed%sUnidentifiedPion",cutName->Data()),"Number of Reconstructed unidentified pions particles");
1873       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedPionAssumingPion",cutName->Data()),"Reconstructed E_{T} from unidentified pions particles assuming pion mass");
1874       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedPionAssumingKaon",cutName->Data()),"Reconstructed E_{T} from unidentified pions particles assuming pion mass");
1875       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedPionAssumingProton",cutName->Data()),"Reconstructed E_{T} from unidentified pions particles assuming pion mass");
1876       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedPion",cutName->Data()),"Reconstructed E_{T} from unidentified pions particles assuming pion mass");
1877
1878       CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentified",cutName->Data()),"Reconstructed E_{T} from unidentified particles using real mass");
1879       CreateEtaPtHisto2D(Form("EtReconstructed%sMisidentifiedElectrons",cutName->Data()),"Reconstructed E_{T} from misidentified electrons");
1880
1881
1882       CreateEtaPtHisto2D(Form("EtReconstructed%sPiPlus",cutName->Data()),"Reconstructed E_{T} from #pi^{+}");
1883       CreateEtaPtHisto2D(Form("EtReconstructed%sPiMinus",cutName->Data()),"Reconstructed E_{T} from #pi^{-}");
1884       CreateEtaPtHisto2D(Form("EtReconstructed%sKPlus",cutName->Data()),"Reconstructed E_{T} from K^{+}");
1885       CreateEtaPtHisto2D(Form("EtReconstructed%sKMinus",cutName->Data()),"Reconstructed E_{T} from K^{-}");
1886       CreateEtaPtHisto2D(Form("EtReconstructed%sProton",cutName->Data()),"Reconstructed E_{T} from p");
1887       CreateEtaPtHisto2D(Form("EtReconstructed%sAntiProton",cutName->Data()),"Reconstructed E_{T} from #bar{p}");
1888       if(fBaryonEnhancement){
1889         CreateEtaPtHisto2D(Form("EtReconstructed%sProtonEnhanced",cutName->Data()),"Reconstructed E_{T} from p");
1890         CreateEtaPtHisto2D(Form("EtReconstructed%sAntiProtonEnhanced",cutName->Data()),"Reconstructed E_{T} from #bar{p}");
1891       }
1892       CreateEtaPtHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),"Reconstructed E_{T} from charged hadrons");
1893       CreateEtaPtHisto2D(Form("EtNReconstructed%sPiPlus",cutName->Data()),"Reconstructed E_{T} from #pi^{+}");
1894       CreateEtaPtHisto2D(Form("EtNReconstructed%sPiMinus",cutName->Data()),"Reconstructed E_{T} from #pi^{-}");
1895       CreateEtaPtHisto2D(Form("EtNReconstructed%sKPlus",cutName->Data()),"Reconstructed E_{T} from K^{+}");
1896       CreateEtaPtHisto2D(Form("EtNReconstructed%sKMinus",cutName->Data()),"Reconstructed E_{T} from K^{-}");
1897       CreateEtaPtHisto2D(Form("EtNReconstructed%sProton",cutName->Data()),"Reconstructed E_{T} from p");
1898       CreateEtaPtHisto2D(Form("EtNReconstructed%sAntiProton",cutName->Data()),"Reconstructed E_{T} from #bar{p}");
1899       if(fBaryonEnhancement){
1900         CreateEtaPtHisto2D(Form("EtNReconstructed%sProtonEnhanced",cutName->Data()),"Reconstructed E_{T} from p");
1901         CreateEtaPtHisto2D(Form("EtNReconstructed%sAntiProtonEnhanced",cutName->Data()),"Reconstructed E_{T} from #bar{p}");
1902       }
1903       CreateEtaPtHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),"Reconstructed E_{T} from charged hadrons");
1904       if(fDataSet==20100){//If this is Pb+Pb
1905         Int_t width = 5;
1906         if(fNCentBins<21) width = 10;
1907         for(Int_t j=0;j<fNCentBins;j++){
1908           CreateEtaPtHisto2D(Form("EtNReconstructed%sPiPlusCB%i",cutName->Data(),j),Form("Reconstructed E_{T} from #pi^{+} for %i-%i central",j*width,(j+1)*width));
1909           CreateEtaPtHisto2D(Form("EtNReconstructed%sPiMinusCB%i",cutName->Data(),j),Form("Reconstructed E_{T} from #pi^{-} for %i-%i central",j*width,(j+1)*width));
1910           CreateEtaPtHisto2D(Form("EtNReconstructed%sKPlusCB%i",cutName->Data(),j),Form("Reconstructed E_{T} from K^{+} for %i-%i central",j*width,(j+1)*width));
1911           CreateEtaPtHisto2D(Form("EtNReconstructed%sKMinusCB%i",cutName->Data(),j),Form("Reconstructed E_{T} from K^{-} for %i-%i central",j*width,(j+1)*width));
1912           CreateEtaPtHisto2D(Form("EtNReconstructed%sProtonCB%i",cutName->Data(),j),Form("Reconstructed E_{T} from p for %i-%i central",j*width,(j+1)*width));
1913           CreateEtaPtHisto2D(Form("EtNReconstructed%sAntiProtonCB%i",cutName->Data(),j),Form("Reconstructed E_{T} from #bar{p} for %i-%i central",j*width,(j+1)*width));
1914           CreateEtaPtHisto2D(Form("EtNReconstructed%sChargedHadronCB%i",cutName->Data(),j),Form("Reconstructed E_{T} from charged hadrons for %i-%i central",j*width,(j+1)*width));
1915         }
1916       }
1917
1918       CreateEtaPtHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),"Reconstructed E_{T} from charged hadrons assuming they are all pions");
1919       CreateEtaPtHisto2D(Form("EtReconstructed%sChargedHadronAssumingProton",cutName->Data()),"Reconstructed E_{T} from charged hadrons assuming they are all pions");
1920       CreateEtaPtHisto2D(Form("EtReconstructed%sChargedHadronAssumingKaon",cutName->Data()),"Reconstructed E_{T} from charged hadrons assuming they are all pions");
1921       CreateEtaPtHisto2D(Form("EtReconstructed%sKPlusAssumingPion",cutName->Data()),"Reconstructed E_{T} from K^{+} assuming #pi mass");
1922       CreateEtaPtHisto2D(Form("EtReconstructed%sKMinusAssumingPion",cutName->Data()),"Reconstructed E_{T} from K^{-} assuming #pi mass");
1923       CreateEtaPtHisto2D(Form("EtReconstructed%sProtonAssumingPion",cutName->Data()),"Reconstructed E_{T} from p assuming #pi mass");
1924       CreateEtaPtHisto2D(Form("EtReconstructed%sAntiProtonAssumingPion",cutName->Data()),"Reconstructed E_{T} from #bar{p} assuming #pi mass");
1925       CreateEtaPtHisto2D(Form("EtReconstructed%sPiPlusAssumingKaon",cutName->Data()),"Reconstructed E_{T} from K^{+} assuming #pi mass");
1926       CreateEtaPtHisto2D(Form("EtReconstructed%sPiMinusAssumingKaon",cutName->Data()),"Reconstructed E_{T} from K^{-} assuming #pi mass");
1927       CreateEtaPtHisto2D(Form("EtReconstructed%sKPlusAssumingKaon",cutName->Data()),"Reconstructed E_{T} from K^{+} assuming #pi mass");
1928       CreateEtaPtHisto2D(Form("EtReconstructed%sKMinusAssumingKaon",cutName->Data()),"Reconstructed E_{T} from K^{-} assuming #pi mass");
1929       CreateEtaPtHisto2D(Form("EtReconstructed%sProtonAssumingKaon",cutName->Data()),"Reconstructed E_{T} from p assuming #pi mass");
1930       CreateEtaPtHisto2D(Form("EtReconstructed%sAntiProtonAssumingKaon",cutName->Data()),"Reconstructed E_{T} from #bar{p} assuming #pi mass");
1931
1932       CreateEtaPtHisto2D(Form("EtReconstructed%sKPlusAssumingProton",cutName->Data()),"Reconstructed E_{T} from K^{+} assuming #pi mass");
1933       CreateEtaPtHisto2D(Form("EtReconstructed%sKMinusAssumingProton",cutName->Data()),"Reconstructed E_{T} from K^{-} assuming #pi mass");
1934       CreateEtaPtHisto2D(Form("EtReconstructed%sPiMinusAssumingProton",cutName->Data()),"Reconstructed E_{T} from p assuming #pi mass");
1935       CreateEtaPtHisto2D(Form("EtReconstructed%sPiPlusAssumingProton",cutName->Data()),"Reconstructed E_{T} from #bar{p} assuming #pi mass");
1936       CreateEtaPtHisto2D(Form("EtReconstructed%sProtonAssumingProton",cutName->Data()),"Reconstructed E_{T} from p assuming #pi mass");
1937       CreateEtaPtHisto2D(Form("EtReconstructed%sAntiProtonAssumingProton",cutName->Data()),"Reconstructed E_{T} from #bar{p} assuming #pi mass");
1938
1939       if(fBaryonEnhancement){
1940         CreateEtaPtHisto2D(Form("EtReconstructed%sProtonAssumingPionEnhanced",cutName->Data()),"Reconstructed E_{T} from p assuming #pi mass");
1941         CreateEtaPtHisto2D(Form("EtReconstructed%sAntiProtonAssumingPionEnhanced",cutName->Data()),"Reconstructed E_{T} from #bar{p} assuming #pi mass");
1942       }
1943
1944       CreateEtaPtHisto2D(Form("EtReconstructed%sEPlus",cutName->Data()),"Reconstructed E_{T} from e^{+}");
1945       CreateEtaPtHisto2D(Form("EtReconstructed%sEMinus",cutName->Data()),"Reconstructed E_{T} from e^{-}");
1946
1947
1948
1949       CreateEtaPtHisto2D(Form("EtReconstructed%sLambdaDaughters",cutName->Data()),"Reconstructed E_{T} from #Lambda Daughters");
1950       CreateEtaPtHisto2D(Form("EtReconstructed%sAntiLambdaDaughters",cutName->Data()),"Reconstructed E_{T} from #bar{#Lambda} Daughters");
1951       CreateEtaPtHisto2D(Form("EtReconstructed%sK0SDaughters",cutName->Data()),"Reconstructed E_{T} from K^{0}_{S} Daughters");
1952       CreateEtaPtHisto2D(Form("EtReconstructed%sLambdaDaughtersReweighted",cutName->Data()),"Reconstructed E_{T} from #Lambda Daughters");
1953       CreateEtaPtHisto2D(Form("EtReconstructed%sAntiLambdaDaughtersReweighted",cutName->Data()),"Reconstructed E_{T} from #bar{#Lambda} Daughters");
1954       CreateEtaPtHisto2D(Form("EtReconstructed%sK0SDaughtersReweighted",cutName->Data()),"Reconstructed E_{T} from K^{0}_{S} Daughters");
1955       CreateEtaPtHisto2D(Form("EtReconstructed%sOmegaDaughters",cutName->Data()),"Reconstructed E_{T} from #Omega^{-} Daughters");
1956       CreateEtaPtHisto2D(Form("EtReconstructed%sAntiOmegaDaughters",cutName->Data()),"Reconstructed E_{T} from #Omega^{+} Daughters");
1957       CreateEtaPtHisto2D(Form("EtReconstructed%sXiDaughters",cutName->Data()),"Reconstructed E_{T} from #Xi^{-} Daughters");
1958       CreateEtaPtHisto2D(Form("EtReconstructed%sAntiXiDaughters",cutName->Data()),"Reconstructed E_{T} from #Xi^{+} Daughters");
1959       CreateEtaPtHisto2D(Form("EtReconstructed%sConversionElectrons",cutName->Data()),"Reconstructed E_{T} from conversion electrons");
1960       CreateEtaPtHisto2D(Form("EtReconstructed%sSecondaryMuons",cutName->Data()),"Reconstructed E_{T} from secondary muons");//from pions
1961       CreateEtaPtHisto2D(Form("EtReconstructed%sSecondaryPions",cutName->Data()),"Reconstructed E_{T} from secondary pions");//from rescattering and sigma+-
1962       CreateEtaPtHisto2D(Form("EtReconstructed%sSecondaryProtons",cutName->Data()),"Reconstructed E_{T} from secondary protons");//from rescattering and sigma+-
1963     }
1964     CreateIntHisto1D(Form("UnidentifiedPIDs%s",cutName->Data()),"PIDs of unidentified particles", "PID", "Number of particles",9, -4,4);
1965     CreateHisto2D(Form("MisidentifiedPIDs%s",cutName->Data()),"PIDs of misidentified particles", "PID real","PID identified",5, -.5,4.5,5, -.5,4.5);
1966     CreateHisto2D(Form("dEdxAll%s",cutName->Data()),"dE/dx for all particles","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
1967     CreateHisto2D(Form("dEdxPion%s",cutName->Data()),"dE/dx for #pi^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
1968     CreateHisto2D(Form("dEdxKaon%s",cutName->Data()),"dE/dx for K^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
1969     CreateHisto2D(Form("dEdxProton%s",cutName->Data()),"dE/dx for p(#bar{p})","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
1970     CreateHisto2D(Form("dEdxElectron%s",cutName->Data()),"dE/dx for e^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
1971     CreateHisto2D(Form("dEdxUnidentified%s",cutName->Data()),"dE/dx for unidentified particles","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
1972   }
1973   delete strTPC;
1974   delete strITS;
1975   delete strTPCITS;
1976
1977   Float_t minEt = 0.0;
1978   Float_t maxEt = 100.0;
1979   Float_t minEtPiKP = 0.0;
1980   Float_t maxEtPiKP = 50.0;
1981   if(fDataSet==20100) maxEt=4000.0;
1982   Int_t nbinsEt = 100;
1983   char histoname[200];
1984   char histotitle[200];
1985   char xtitle[50];
1986   char ytitle[50];
1987   TString *sTPC = new TString("TPC");
1988   TString *sITS = new TString("ITS");
1989   TString *sTPCpt = new TString("0.15");
1990   TString *sITSpt = new TString("0.10");
1991   TString *sPID = new TString("");
1992   TString *sNoPID = new TString("NoPID");
1993   TString *sNoPIDString = new TString(", No PID");
1994   TString *sHadEt = new TString("HadEt");
1995   TString *sTotEt = new TString("TotEt");
1996   TString *sPiKPEt =  new TString("PiKPEt");
1997   TString *sTotEtString = new TString("total E_{T}");
1998   TString *sHadEtString = new TString("hadronic E_{T}");
1999   TString *sPiKPEtString = new TString("E_{T}^{#pi,K,p}");
2000   TString *sFull = new TString("Full");
2001   TString *sEMCAL = new TString("EMCAL");
2002   TString *sPHOS = new TString("PHOS");
2003   float etDiff = 1.5;
2004   float etDiffLow = etDiff;
2005   if(fDataSet!=20100){//If this is p+p
2006     etDiffLow = 2.5;
2007   }
2008
2009   for(int tpc = 0;tpc<lastcutset;tpc++){
2010     TString *detector = NULL;
2011     TString *ptstring = NULL;
2012     if(tpc==1) {detector = sTPC; ptstring = sTPCpt;}
2013     else{detector = sITS; ptstring = sITSpt;}
2014     for(int hadet = 0;hadet<2;hadet++){
2015       TString *et = NULL;
2016       TString *etstring = NULL;
2017       if(hadet==1) {et = sHadEt; etstring = sHadEtString;}
2018       else{et = sTotEt; etstring = sTotEtString;}
2019       for(int type = 0;type<3;type++){
2020         if(type==0 && !fInvestigateFull) continue;
2021         if(type==1 && !fInvestigateEMCal) continue;
2022         if(type==2 && !fInvestigatePHOS) continue;
2023         TString *acceptance = NULL;
2024         switch(type){
2025         case 0:
2026           acceptance = sFull;
2027           etDiff = 1.5;
2028           break;
2029         case 1:
2030           acceptance = sEMCAL;
2031           etDiff = 5;
2032           break;
2033         case 2:
2034           acceptance = sPHOS;
2035           etDiff = 5;
2036           break;
2037         default:
2038           acceptance = sFull;
2039         }
2040         for(int pid = 0;pid<2;pid++){
2041           TString *partid = NULL;
2042           TString *partidstring = NULL;
2043           if(pid==1){partid = sPID; partidstring = sPID;}
2044           else{partid = sNoPID; partidstring = sNoPIDString;}
2045
2046           snprintf(histoname,200,"Sim%sMinusReco%s%sAcceptance%s%s",et->Data(),et->Data(),acceptance->Data(),detector->Data(),partid->Data());
2047           snprintf(histotitle,200,"(Simulated %s - reconstructed %s)/(Simulated %s) with %s acceptance for p_{T}>%s GeV/c%s",etstring->Data(),etstring->Data(),etstring->Data(),acceptance->Data(),ptstring->Data(),partidstring->Data());
2048           snprintf(ytitle,50,"(Simulated %s - reconstructed %s)/(Simulated %s)",etstring->Data(),etstring->Data(),etstring->Data());
2049           snprintf(xtitle,50,"Simulated %s",etstring->Data());
2050           CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiffLow,etDiff);
2051
2052           if(type==0){//only fill for full acceptance
2053             snprintf(histoname,200,"Sim%sVsReco%s%sAcceptance%s%s",et->Data(),et->Data(),acceptance->Data(),detector->Data(),partid->Data());
2054             snprintf(histotitle,200,"Simulated %s vs reconstructed %s with %s acceptance for p_{T}>%s GeV/c%s",etstring->Data(),etstring->Data(),acceptance->Data(),ptstring->Data(),partidstring->Data());
2055             snprintf(ytitle,50,"Reconstructed %s",etstring->Data());
2056             snprintf(xtitle,50,"Simulated %s",etstring->Data());
2057             CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt*4,minEt,maxEt,nbinsEt*4,minEt,maxEt);
2058             if(hadet==1){//on top of that we want to create pikp histograms without adding a full option parallel to had and tot et - therefore we will just create these histos when we create the hadet histos
2059               snprintf(histoname,200,"Sim%sVsReco%s%sAcceptance%s%s",sPiKPEt->Data(),sPiKPEt->Data(),acceptance->Data(),detector->Data(),partid->Data());
2060               snprintf(histotitle,200,"Simulated %s vs reconstructed %s with %s acceptance for p_{T}>%s GeV/c%s",sPiKPEtString->Data(),sPiKPEtString->Data(),acceptance->Data(),ptstring->Data(),partidstring->Data());
2061               snprintf(ytitle,50,"Reconstructed %s",sPiKPEtString->Data());
2062               snprintf(xtitle,50,"Simulated %s",sPiKPEtString->Data());
2063               CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt*4,minEtPiKP,maxEtPiKP,nbinsEt*4,minEtPiKP,maxEtPiKP);
2064             }
2065           }
2066
2067           if(hadet==0 && type==0 && fInvestigatePiKP){//we only want to do this once...  not the most elegant way of coding but hey...
2068             snprintf(histoname,200,"SimPiKPMinusRecoPiKP%sAcceptance%s%s",acceptance->Data(),detector->Data(),partid->Data());
2069             snprintf(histotitle,200,"(Sim PiKP - reco PiKP)/(Sim PiKP) with %s acceptance for p_{T}>%s GeV/c%s",acceptance->Data(),ptstring->Data(),partidstring->Data());
2070             snprintf(ytitle,50,"(Sim PiKP - reco PiKP)/(Sim PiKP)");
2071             snprintf(xtitle,50,"Simulated E_{T}^{#pi,K,p}");
2072             CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiffLow,etDiff);
2073           }
2074         }
2075       }
2076     }
2077   }
2078   CreateHisto1D("SimPiKPEt","Simulated #pi,K,p E_{T}","Simulated #pi,K,p E_{T}","Number of events",nbinsEt*4,minEtPiKP,maxEtPiKP);
2079   CreateHisto1D("SimPiKPEtND","Simulated #pi,K,p E_{T}","Simulated #pi,K,p E_{T}","Number of events",nbinsEt*4,minEtPiKP,maxEtPiKP);
2080   CreateHisto1D("SimPiKPEtDD","Simulated #pi,K,p E_{T}","Simulated #pi,K,p E_{T}","Number of events",nbinsEt*4,minEtPiKP,maxEtPiKP);
2081   CreateHisto1D("SimPiKPEtSD","Simulated #pi,K,p E_{T}","Simulated #pi,K,p E_{T}","Number of events",nbinsEt*4,minEtPiKP,maxEtPiKP);
2082   CreateHisto1D("SimTotEt","Simulated Total E_{T}","Simulated Total E_{T}","Number of events",nbinsEt*4,minEt,maxEt);
2083   CreateHisto1D("SimHadEt","Simulated Hadronic E_{T}","Simulated Hadronic E_{T}","Number of events",nbinsEt*4,minEt,maxEt);
2084   CreateHisto1D("SimTotEtND","Simulated Total E_{T}","Simulated Total E_{T} for non-diffractive events","Number of events",nbinsEt*4,minEt,maxEt);
2085   CreateHisto1D("SimHadEtND","Simulated Hadronic E_{T}","Simulated Hadronic E_{T} for non-diffractive events","Number of events",nbinsEt*4,minEt,maxEt);
2086   CreateHisto1D("SimTotEtSD","Simulated Total E_{T}","Simulated Total E_{T} for singly diffractive events","Number of events",nbinsEt*4,minEt,maxEt);
2087   CreateHisto1D("SimHadEtSD","Simulated Hadronic E_{T}","Simulated Hadronic E_{T} for singly diffractive events","Number of events",nbinsEt*4,minEt,maxEt);
2088   CreateHisto1D("SimTotEtDD","Simulated Total E_{T}","Simulated Total E_{T} for doubly diffractive events","Number of events",nbinsEt*4,minEt,maxEt);
2089   CreateHisto1D("SimHadEtDD","Simulated Hadronic E_{T}","Simulated Hadronic E_{T} for doubly diffractive events","Number of events",nbinsEt*4,minEt,maxEt);
2090   if(fDataSet==20100){
2091     Int_t width = 5;
2092     if(fNCentBins<21) width = 10;
2093     for(Int_t j=0;j<fNCentBins;j++){
2094       CreateHisto1D(Form("SimTotEtCB%i",j),Form("Simulated Total E_{T} for %i-%i central",j*width,(j+1)*width),"Simulated Total E_{T}","Number of events",nbinsEt*4,minEt,maxEt);
2095       CreateHisto1D(Form("SimHadEtCB%i",j),Form("Simulated Hadronic E_{T} for %i-%i central",j*width,(j+1)*width),"Simulated Hadronic E_{T}","Number of events",nbinsEt*4,minEt,maxEt);
2096       CreateHisto1D(Form("SimPiKPEtCB%i",j),Form("Simulated #pi,K,p E_{T} for %i-%i central",j*width,(j+1)*width),"Simulated #pi,K,p E_{T}","Number of events",nbinsEt,minEt,maxEt);
2097     }
2098   }
2099
2100   etDiff = 0.15;
2101
2102   if(fInvestigateSmearing && !fRunLightweight){
2103     //======================================================================
2104
2105     snprintf(histoname,200,"SimPiKPEtMeasMinusEtRealPiKP");
2106     snprintf(histotitle,200,"Simulated (all reconstructed - primaries)/all reconstructed for reconstructed tracks only");
2107     snprintf(ytitle,50,"(primary-all)/primary");
2108     snprintf(xtitle,50,"true p, K, p E_{T} for primary tracks");
2109     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff*5,etDiff*5);
2110
2111     snprintf(histoname,200,"SimPiKPEtMinusSimAllCorrSmearedRecoOnly");
2112     snprintf(histotitle,200,"Simulated (primary-all)/primary for reconstructed tracks only");
2113     snprintf(ytitle,50,"(primary-all)/primary");
2114     snprintf(xtitle,50,"true p, K, p E_{T} for primary tracks");
2115     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff*5,etDiff*5);
2116
2117     snprintf(histoname,200,"SimPiKPEtMinusSimEffCorrRecoOnly");
2118     snprintf(histotitle,200,"(sim-reco)/sim primary #pi,k,p for p_{T}>0.15");
2119     snprintf(ytitle,50,"(sim-reco)/sim");
2120     snprintf(xtitle,50,"true p, K, p E_{T} for primary tracks");
2121     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff*5,etDiff*5);
2122
2123     snprintf(histoname,200,"SimPiKPEtMinusSimEffBkgdCorrRecoOnly");
2124     snprintf(histotitle,200,"(sim-reco)/sim primary #pi,k,p for p_{T}>0.15 with background subtraction");
2125     snprintf(ytitle,50,"(sim-reco)/sim");
2126     snprintf(xtitle,50,"true p, K, p E_{T} for primary tracks");
2127     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff*5,etDiff*5);
2128
2129     snprintf(histoname,200,"SimPiKPEtMinusSimEffCorrRecoPiOnly");
2130     snprintf(histotitle,200,"(sim-reco)/sim primary #pi for p_{T}>0.15");
2131     snprintf(ytitle,50,"(sim-reco)/sim");
2132     snprintf(xtitle,50,"true #pi E_{T} for primary tracks");
2133     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt/2,nbinsEt,-etDiff*5,etDiff*5);
2134
2135     snprintf(histoname,200,"SimPiKPEtMinusSimEffCorrRecoKOnly");
2136     snprintf(histotitle,200,"(sim-reco)/sim primary K for p_{T}>0.15");
2137     snprintf(ytitle,50,"(sim-reco)/sim");
2138     snprintf(xtitle,50,"true K E_{T} for primary tracks");
2139     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt/6,nbinsEt,-etDiff*5,etDiff*5);
2140
2141     snprintf(histoname,200,"SimPiKPEtMinusSimEffCorrRecoPOnly");
2142     snprintf(histotitle,200,"(sim-reco)/sim primary p for p_{T}>0.15");
2143     snprintf(ytitle,50,"(sim-reco)/sim");
2144     snprintf(xtitle,50,"true p E_{T} for primary tracks");
2145     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt/6,nbinsEt,-etDiff*5,etDiff*5);
2146
2147     snprintf(histoname,200,"SimPiKPEtMinusSimAllSmearedRecoOnly");
2148     snprintf(histotitle,200,"Simulated (primary-all)/primary for reconstructed tracks only");
2149     snprintf(ytitle,50,"(primary-all)/primary");
2150     snprintf(xtitle,50,"true p, K, p E_{T} for primary tracks");
2151     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff*5,0.0);
2152
2153     snprintf(histoname,200,"SimPiKPEtMinusSimPIDSmearedRecoOnly");
2154     snprintf(histotitle,200,"Simulated (true-smeared)/true for reconstructed tracks only with PID smearing");
2155     snprintf(ytitle,50,"(true-smeared)/true");
2156     snprintf(xtitle,50,"true p, K, p E_{T}");
2157     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff,etDiff);
2158
2159     snprintf(histoname,200,"SimPiKPEtMinusSimSmearedRecoOnly");
2160     snprintf(histotitle,200,"Simulated (true-smeared)/true for reconstructed tracks only");
2161     snprintf(ytitle,50,"(true-smeared)/true");
2162     snprintf(xtitle,50,"true p, K, p E_{T}");
2163     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff/15,etDiff/15);
2164
2165     snprintf(histoname,200,"SimPiKPPtMinusSimSmearedRecoOnly");
2166     snprintf(histotitle,200,"Simulated (true-smeared)/true for reconstructed tracks only");
2167     snprintf(ytitle,50,"(true-smeared)/true");
2168     snprintf(xtitle,50,"true p, K, p p_{T}");
2169     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff/15,etDiff/15);
2170
2171     snprintf(histoname,200,"SimPiKPEtMinusSimSmearedMultRecoOnly");
2172     snprintf(histotitle,200,"Simulated (true-smeared)/true for reconstructed tracks only");
2173     snprintf(ytitle,50,"(true-smeared)/true");
2174     snprintf(xtitle,50,"number of reconstructed particles");
2175     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff/15,etDiff/15);
2176
2177     snprintf(histoname,200,"SimPiKPPtMinusSimSmearedMultRecoOnly");
2178     snprintf(histotitle,200,"Simulated (true-smeared)/true for reconstructed tracks only");
2179     snprintf(ytitle,50,"(true-smeared)/true");
2180     snprintf(xtitle,50,"number of reconstructed particles");
2181     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff/15,etDiff/15);
2182
2183     //======================================================================
2184
2185     snprintf(histoname,200,"SimPiKPEtMinusSimPtSmeared");
2186     snprintf(histotitle,200,"Simulated (true-smeared)/true for 0.5 percent momentum smearing");
2187     snprintf(ytitle,50,"(true-smeared)/true");
2188     snprintf(xtitle,50,"true p, K, p E_{T}");
2189     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff,etDiff);
2190     snprintf(histoname,200,"SimPiKPEtPtSmeared");
2191     snprintf(histotitle,200,"Simulated E_{T} for 0.5 percent momentum smearing");
2192     snprintf(ytitle,50,"Number of events");
2193     snprintf(xtitle,50,"p, K, p E_{T}");
2194     CreateHisto1D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt);
2195
2196     //======================================================================
2197
2198     snprintf(histoname,200,"SimPiKPEtMinusSimEfficiencySmeared");
2199     snprintf(histotitle,200,"Simulated (true-smeared)/true for efficiency smearing");
2200     snprintf(ytitle,50,"(true-smeared)/true");
2201     snprintf(xtitle,50,"true p, K, p E_{T}");
2202     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff*7,etDiff*7);
2203     snprintf(histoname,200,"SimPiKPEtEfficiencySmeared");
2204     snprintf(histotitle,200,"Simulated E_{T} for efficiency smearing");
2205     snprintf(ytitle,50,"Number of events");
2206     snprintf(xtitle,50,"p, K, p E_{T}");
2207     CreateHisto1D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt);
2208
2209     //======================================================================
2210
2211     snprintf(histoname,200,"SimPiKPEtMinusSimPtCutSmearedTPC");
2212     snprintf(histotitle,200,"Simulated (true-smeared)/true for p_{T}>0.15 GeV/c smearing");
2213     snprintf(ytitle,50,"(true-smeared)/true");
2214     snprintf(xtitle,50,"true p, K, p E_{T}");
2215     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff,etDiff);
2216     snprintf(histoname,200,"SimPiKPEtPtCutSmearedTPC");
2217     snprintf(histotitle,200,"Simulated E_{T} for p_{T}>0.15 GeV/c smearing");
2218     snprintf(ytitle,50,"Number of events");
2219     snprintf(xtitle,50,"p, K, p E_{T}");
2220     CreateHisto1D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt);
2221
2222
2223     //======================================================================
2224
2225     snprintf(histoname,200,"SimPiKPEtMinusSimPtCutSmearedITS");
2226     snprintf(histotitle,200,"Simulated (true-smeared)/true for p_{T}>0.10 GeV/c smearing");
2227     snprintf(ytitle,50,"(true-smeared)/true");
2228     snprintf(xtitle,50,"true p, K, p E_{T}");
2229     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff,etDiff);
2230     snprintf(histoname,200,"SimPiKPEtPtCutSmearedITS");
2231     snprintf(histotitle,200,"Simulated E_{T} for p_{T}>0.10 GeV/c smearing");
2232     snprintf(ytitle,50,"Number of events");
2233     snprintf(xtitle,50,"p, K, p E_{T}");
2234     CreateHisto1D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt);
2235
2236     //======================================================================
2237
2238     snprintf(histoname,200,"SimPiKPEtMinusSimPIDSmeared");
2239     snprintf(histotitle,200,"Simulated (true-smeared)/true for PID smearing");
2240     snprintf(ytitle,50,"(true-smeared)/true");
2241     snprintf(xtitle,50,"true p, K, p E_{T}");
2242     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff,etDiff);
2243     snprintf(histoname,200,"SimPiKPEtPIDSmeared");
2244     snprintf(histotitle,200,"Simulated E_{T} for PID smearing");
2245     snprintf(ytitle,50,"Number of events");
2246     snprintf(xtitle,50,"p, K, p E_{T}");
2247     CreateHisto1D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt);
2248
2249     //======================================================================
2250
2251     snprintf(histoname,200,"SimPiKPEtMinusSimPIDSmearedNoID");
2252     snprintf(histotitle,200,"Simulated (true-smeared)/true for PID smearing No ID");
2253     snprintf(ytitle,50,"(true-smeared)/true");
2254     snprintf(xtitle,50,"true p, K, p E_{T}");
2255     CreateHisto2D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt,nbinsEt,-etDiff,etDiff);
2256     snprintf(histoname,200,"SimPiKPEtPIDSmearedNoID");
2257     snprintf(histotitle,200,"Simulated E_{T} for PID smearing No ID");
2258     snprintf(ytitle,50,"Number of events");
2259     snprintf(xtitle,50,"p, K, p E_{T}");
2260     CreateHisto1D(histoname,histotitle,xtitle,ytitle,nbinsEt,minEt,maxEt);
2261   }
2262   delete sTPC;
2263   delete sITS;
2264   delete sTPCpt;
2265   delete sITSpt;
2266   delete sPID;
2267   delete sNoPID;
2268   delete sNoPIDString;
2269   delete sHadEt;
2270   delete sTotEt;
2271   delete sPiKPEt;
2272   delete sTotEtString;
2273   delete sHadEtString;
2274   delete sPiKPEtString;
2275   delete sFull;
2276   delete sEMCAL;
2277   delete sPHOS;
2278   CreateIntHisto1D("NEvents","Number of events","number of events","Number of events",1,0,1);
2279   CreateIntHisto1D("NEventsSD","Number of events","number of singly diffractive events","Number of events",1,0,1);
2280   CreateIntHisto1D("NEventsDD","Number of events","number of doubly diffractive events","Number of events",1,0,1);
2281   CreateIntHisto1D("NEventsND","Number of events","number of non-diffractive events","Number of events",1,0,1);
2282   if( !fRunLightweight){
2283     CreateResolutionPtHisto2D("presolutionTPC","p resolution","p^{rec}","(p^{sim}-p^{rec})/p^{rec}");
2284     CreateResolutionPtHisto2D("pTresolutionTPC","p_{T} resolution","p_{T}^{rec}","(p_{T}^{sim}-p_{T}^{rec})/p_{T}^{rec}");
2285     CreateResolutionPtHisto2D("ETresolutionTPC","E_{T} resolution","E_{T}^{rec}","(E_{T}^{sim}-E_{T}^{rec})/E_{T}^{rec}");
2286     CreateResolutionPtHisto2D("pTresolutionTPCITS","p_{T} resolution","p_{T}^{rec}","(p_{T}^{sim}-p_{T}^{rec})/p_{T}^{rec}");
2287     CreateResolutionPtHisto2D("ETresolutionTPCITS","E_{T} resolution","E_{T}^{rec}","(E_{T}^{sim}-E_{T}^{rec})/E_{T}^{rec}");
2288     CreateResolutionPtHisto2D("presolutionTPCITS","p resolution","p^{rec}","(p^{sim}-p^{rec})/p^{rec}");
2289     CreateResolutionPtHisto2D("pTresolutionITS","p_{T} resolution","p_{T}^{rec}","(p_{T}^{sim}-p_{T}^{rec})/p_{T}^{rec}");
2290     CreateResolutionPtHisto2D("ETresolutionITS","E_{T} resolution","E_{T}^{rec}","(E_{T}^{sim}-E_{T}^{rec})/E_{T}^{rec}");
2291     CreateResolutionPtHisto2D("presolutionITS","p resolution","p^{rec}","(p^{sim}-p^{rec})/p^{rec}");
2292     CreatePtHisto1D("pTsimITS","p_{T}^{sim}","p_{T}^{sim}","Number of particles");
2293     CreatePtHisto1D("pTsimTPC","p_{T}^{sim}","p_{T}^{sim}","Number of particles");
2294     CreatePtHisto1D("pTsimTPCITS","p_{T}^{sim}","p_{T}^{sim}","Number of particles");
2295     CreatePtHisto1D("pTrecITS","p_{T}^{rec}","p_{T}^{rec}","Number of particles");
2296     CreatePtHisto1D("pTrecTPC","p_{T}^{rec}","p_{T}^{rec}","Number of particles");
2297     CreatePtHisto1D("pTrecTPCITS","p_{T}^{rec}","p_{T}^{rec}","Number of particles");
2298     if(fDataSet==20100){
2299       Int_t width = 5;
2300       if(fNCentBins<21) width = 10;
2301       for(Int_t j=0;j<fNCentBins;j++){
2302         CreatePtHisto1D(Form("pTsimITSCB%i",j),Form("p_{T}^{sim} for %i-%i central",j*width,(j+1)*width),"p_{T}^{sim}","Number of particles");
2303         CreatePtHisto1D(Form("pTsimTPCITSCB%i",j),Form("p_{T}^{sim} for %i-%i central",j*width,(j+1)*width),"p_{T}^{sim}","Number of particles");
2304         CreatePtHisto1D(Form("pTsimTPCCB%i",j),Form("p_{T}^{sim} for %i-%i central",j*width,(j+1)*width),"p_{T}^{sim}","Number of particles");
2305         CreatePtHisto1D(Form("pTrecITSCB%i",j),Form("p_{T}^{rec} for %i-%i central",j*width,(j+1)*width),"p_{T}^{rec}","Number of particles");
2306         CreatePtHisto1D(Form("pTrecTPCITSCB%i",j),Form("p_{T}^{rec} for %i-%i central",j*width,(j+1)*width),"p_{T}^{rec}","Number of particles");
2307         CreatePtHisto1D(Form("pTrecTPCCB%i",j),Form("p_{T}^{rec} for %i-%i central",j*width,(j+1)*width),"p_{T}^{rec}","Number of particles");
2308       }
2309     }
2310   }
2311
2312 }
2313