PWGLF/totEt/AliAnalysisEt.cxx

   1 //_________________________________________________________________________
   2 //  Utility Class for transverse energy studies
   3 //  Base class for ESD & MC analysis
   4 //  - reconstruction and MonteCarlo output
   5 // implementation file
   6 //
   7 //*-- Authors: Oystein Djuvsland (Bergen), David Silvermyr (ORNL)
   8 //_________________________________________________________________________
   9
  10 #include "AliAnalysisEt.h"
  11 #include "TMath.h"
  12 #include "TList.h"
  13 #include "TH1F.h"
  14 #include "TH2F.h"
  15 #include "TH1I.h"
  16 #include "TTree.h"
  17 #include <iostream>
  18 #include "AliAnalysisEtCuts.h"
  19 #include "AliESDtrackCuts.h"
  20 #include "AliESDCaloCluster.h"
  21 #include "AliVEvent.h"
  22 #include "Rtypes.h"
  23 #include "TString.h"
  24 #include "AliCentrality.h"
  25 #include "AliAnalysisEtSelector.h"
  26 #include "AliAnalysisEtTrackMatchCorrections.h"
  27 #include "AliAnalysisEtRecEffCorrection.h"
  28 #include "TFile.h"
  29 #include "TVector3.h"
  30 #include "AliPIDResponse.h"
  31 #include "AliTPCPIDResponse.h"
  32 #include "AliInputEventHandler.h"
  33 #include "AliAnalysisManager.h"
  34
  35 using namespace std;
  36 ClassImp(AliAnalysisEt);
  37
  38 /* Auxiliary Histogram variables */
  39 Int_t AliAnalysisEt::fgnumOfEtaBins = 16;
  40 Float_t AliAnalysisEt::fgEtaAxis[17]={-0.78, -0.7, -0.58, -0.46, -0.34, -0.22, -0.12, -0.06, -0.0, 0.06, 0.12, 0.22, 0.34, 0.46, 0.58, 0.7, 0.78};
  41
  42 Int_t AliAnalysisEt::fgNumOfPtBins = 111;
  43 Float_t AliAnalysisEt::fgPtAxis[117]=
  44    {0.0,0.01,0.02,0.03,0.04, 0.05, 0.06,0.07,0.08,0.09, 0.10,0.11, .12,0.13, .14,0.15, .16,0.17, .18,0.19,
  45         0.2, .22, .24, .26, .28, 0.30, 0.32, .34, .36, .38, 0.40, .42, .44, .46, .48,
  46         0.5, .52, .54, .56, .58, 0.60, 0.62, .64, .66, .68, 0.70, .72, .74, .76, .78,
  47         .80, .82, .84, .86, .88, 0.90, 0.92, .94, .96, .98, 1.00,1.05, 1.1,1.15, 1.2,
  48         1.25, 1.3,1.35,1.40,1.45, 1.50, 1.55, 1.6,1.65, 1.7, 1.75, 1.8,1.85, 1.9,1.95,
  49         2.0, 2.2, 2.4, 2.6, 2.8, 3.00, 3.20, 3.4, 3.6, 3.8, 4.00, 4.2, 4.4, 4.6, 4.8,
  50         5.0, 5.5, 6.0, 6.5, 7.0, 7.50, 8.00, 8.5, 9.0, 9.5, 10.0,12.0,14.0,16.0,18.0,
  51         20.0,25.0,30.0,35.0,40.0, 45.0, 50.0};
  52
  53 Int_t AliAnalysisEt::fgNumOfEBins = 78;
  54 Float_t AliAnalysisEt::fgEAxis[79]={0., 0.05, 0.10, 0.15, 0.20, 0.25, 0.30, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95,
  55                                                                         1.,1.5,2.,2.5,3.,3.5,4.,4.5,5.,5.5,6.,6.5,7.,7.5,8.,8.5,9.,9.5,10.,11.,
  56                                                                         12.,13.,14.,15.,16.,17.,18.,19.,20.,21.,22.,23.,24.,25.,26.,27.,28.,29.,30.,31.,
  57                                                                         32.,33.,34.,35.,36.,37.,38.,39.,40.,41.,42.,43.,44.,45.,46.,47.,48.,49.,50.};
  58
  59 Int_t AliAnalysisEt::fgNumOfRBins = 47;
  60 Float_t AliAnalysisEt::fgRAxis[48]={-2.,-1.,0.,0.0005,0.001,0.0015,0.002,0.0025,0.003,0.0035,0.004,0.0045,0.005,0.0055,0.006,0.0065,0.007,0.0075,0.008,0.0085,0.009,0.0095,0.01,
  61                                                                         0.011,0.012,0.013,0.014,0.015,0.016,0.017,0.018,0.019,0.020,0.022,0.024,0.026,0.028,0.03,0.032,0.034,0.036,0.038,0.04,0.042,0.044,0.046,0.048,0.05};
  62
  63
  64 AliAnalysisEt::AliAnalysisEt() : AliAnalysisEtCommon()
  65                                ,fTmCorrections(0)
  66                                ,fReCorrections(0)
  67                                ,fEventSummaryTree(0)
  68                                ,fAcceptedTree(0)
  69                                ,fDepositTree(0)
  70                                ,fTotEt(0)
  71                                ,fTotEtAcc(0)
  72                                ,fTotNeutralEt(0)
  73                                ,fTotNeutralEtAcc(0)
  74                                ,fTotChargedEt(0)
  75                                ,fTotChargedEtAcc(0)
  76                                ,fMultiplicity(0)
  77                                ,fChargedMultiplicity(0)
  78                                ,fNeutralMultiplicity(0)
  79                                ,fProtonEt(0)
  80                                ,fAntiProtonEt(0)
  81                                ,fNeutronEt(0)
  82                                ,fAntiNeutronEt(0)
  83                                ,fPi0Et(0)
  84                                ,fPiPlusEt(0)
  85                                ,fPiMinusEt(0)
  86                                ,fKPlusEt(0)
  87                                ,fKMinusEt(0)
  88                                ,fK0sEt(0)
  89                                ,fK0lEt(0)
  90                                ,fMuMinusEt(0)
  91                                ,fMuPlusEt(0)
  92                                ,fEMinusEt(0)
  93                                ,fEPlusEt(0)
  94                                ,fGammaEt(0)
  95                                ,fProtonRemovedEt(0)
  96                                ,fAntiProtonRemovedEt(0)
  97                                ,fNeutronRemovedEt(0)
  98                                ,fAntiNeutronRemovedEt(0)
  99                                ,fPi0RemovedEt(0)
 100                                ,fPiPlusRemovedEt(0)
 101                                ,fPiMinusRemovedEt(0)
 102                                ,fKPlusRemovedEt(0)
 103                                ,fKMinusRemovedEt(0)
 104                                ,fK0sRemovedEt(0)
 105                                ,fK0lRemovedEt(0)
 106                                ,fMuMinusRemovedEt(0)
 107                                ,fMuPlusRemovedEt(0)
 108                                ,fEMinusRemovedEt(0)
 109                                ,fEPlusRemovedEt(0)
 110                                ,fGammaRemovedEt(0)
 111                                ,fProtonMult(0)
 112                                ,fAntiProtonMult(0)
 113                                ,fNeutronMult(0)
 114                                ,fAntiNeutronMult(0)
 115                                ,fPi0Mult(0)
 116                                ,fPiPlusMult(0)
 117                                ,fPiMinusMult(0)
 118                                ,fKPlusMult(0)
 119                                ,fKMinusMult(0)
 120                                ,fK0sMult(0)
 121                                ,fK0lMult(0)
 122                                ,fMuMinusMult(0)
 123                                ,fMuPlusMult(0)
 124                                ,fEMinusMult(0)
 125                                ,fEPlusMult(0)
 126                                ,fGammaMult(0)
 127                                ,fProtonRemovedMult(0)
 128                                ,fAntiProtonRemovedMult(0)
 129                                ,fNeutronRemovedMult(0)
 130                                ,fAntiNeutronRemovedMult(0)
 131                                ,fPi0RemovedMult(0)
 132                                ,fPiPlusRemovedMult(0)
 133                                ,fPiMinusRemovedMult(0)
 134                                ,fKPlusRemovedMult(0)
 135                                ,fKMinusRemovedMult(0)
 136                                ,fK0sRemovedMult(0)
 137                                ,fK0lRemovedMult(0)
 138                                ,fMuMinusRemovedMult(0)
 139                                ,fMuPlusRemovedMult(0)
 140                                ,fEMinusRemovedMult(0)
 141                                ,fEPlusRemovedMult(0)
 142                                ,fGammaRemovedMult(0)
 143                                ,fEnergyDeposited(0)
 144                                ,fMomentumTPC(0)
 145                                ,fCharge(0)
 146                                ,fParticlePid(0)
 147                                ,fPidProb(0)
 148                                ,fTrackPassedCut(kFALSE)
 149                                ,fCentClass(0)
 150                                ,fDetectorRadius(0)
 151                                ,fSingleCellEnergyCut(0)
 152                                ,fChargedEnergyRemoved(0)
 153                                ,fNeutralEnergyRemoved(0)
 154                                ,fGammaEnergyAdded(0)
 155                                ,fHistEt(0)
 156                                ,fHistNeutralMult(0)
 157                                ,fHistPhivsPtPos(0)
 158                                ,fHistPhivsPtNeg(0)
 159                                ,fCentrality(0)
 160                                ,fMakeSparse(kFALSE)
 161                                ,fCutFlow(0)
 162                                ,fSelector(0)
 163                                ,fPIDResponse(0)
 164                                ,fsub(1.0)
 165                                ,fsubmeanhade(1.0)
 166
 167 {}
 168
 169 AliAnalysisEt::~AliAnalysisEt()
 170 {//Destructor
 171   delete fTmCorrections;
 172   delete fReCorrections;
 173   if(fDepositTree){
 174     fDepositTree->Clear();
 175     delete fDepositTree; // optional TTree
 176   }
 177   if(fEventSummaryTree)
 178   {
 179     fEventSummaryTree->Clear();
 180     delete fEventSummaryTree;
 181   }
 182   if(fAcceptedTree)
 183   {
 184     fAcceptedTree->Clear();
 185     delete fAcceptedTree;
 186   }
 187   delete fHistEt; //Et spectrum
 188   delete fHistNeutralMult; //Neutral multiplicity
 189   delete fHistPhivsPtPos; //phi vs pT plot for positive tracks
 190   delete fHistPhivsPtNeg; //phi vs pT Moplot for negative tracks
 191   //delete fCentrality;//this code does not actually own AliCentrality so we don't have to worry about deleting it...  we just borrow it...
 192   delete fCutFlow;
 193   delete fSelector;
 194   delete fPIDResponse;
 195 }
 196
 197 void AliAnalysisEt::FillOutputList(TList *list)
 198 { // histograms to be added to output
 199     list->Add(fHistEt);
 200     list->Add(fHistNeutralMult);
 201
 202     list->Add(fHistPhivsPtPos);
 203     list->Add(fHistPhivsPtNeg);
 204
 205     if (fCuts) {
 206       if (fCuts->GetHistMakeTree()) {
 207         //list->Add(fTree);
 208         list->Add(fEventSummaryTree);
 209       }
 210       if (fCuts->GetHistMakeTreeDeposit()) {
 211         list->Add(fDepositTree);
 212       }
 213     }
 214
 215     list->Add(fCutFlow);
 216
 217     AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 218     if (!man) {
 219       AliFatal("Analysis manager needed");
 220       return;
 221     }
 222
 223     AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(man->GetInputEventHandler());
 224     if (!inputHandler) {
 225       AliFatal("Input handler needed");
 226       return;
 227     }
 228
 229     //pid response object
 230     fPIDResponse=inputHandler->GetPIDResponse();
 231     if (!fPIDResponse) AliError("PIDResponse object was not created");
 232
 233 }
 234
 235 void AliAnalysisEt::Init()
 236 {// clear variables, set up cuts and PDG info
 237   AliAnalysisEtCommon::Init();
 238   if(ReadCorrections("calocorrections.root") != 0)
 239   {
 240     // Shouldn't do this, why oh why are exceptions not allowed?
 241     exit(-1);
 242   }
 243   ResetEventValues();
 244
 245 }
 246
 247 void AliAnalysisEt::CreateHistograms()
 248 { // create histograms..
 249   // histogram binning for E_T, p_T and Multiplicity: defaults for p+p
 250   Int_t nbinsEt = 10000;
 251   Double_t minEt = 0.0;
 252   Double_t maxEt = 1000;
 253   Int_t nbinsPt = 200;
 254   Double_t minPt = 0;
 255   Double_t maxPt = 20;
 256   Int_t nbinsMult = 200;
 257   Double_t minMult = -0.5; // offset -0.5 to have integer bins centered around 0
 258   Double_t maxMult = nbinsMult + minMult; // 1 bin per integer value
 259
 260   // see if we should change histogram limits etc, and possibly create a tree
 261   if (fCuts) {
 262     if (fCuts->GetHistMakeTree()) {
 263       CreateTrees();
 264     }
 265
 266     nbinsMult = fCuts->GetHistNbinsMult();
 267     minMult = fCuts->GetHistMinMult();
 268     maxMult = fCuts->GetHistMaxMult();
 269
 270     nbinsEt = fCuts->GetHistNbinsTotEt();
 271     minEt = fCuts->GetHistMinTotEt();
 272     maxEt = fCuts->GetHistMaxTotEt();
 273
 274     nbinsPt = fCuts->GetHistNbinsParticlePt();
 275     minPt = fCuts->GetHistMinParticlePt();
 276     maxPt = fCuts->GetHistMaxParticlePt();
 277   }
 278
 279     TString histname = "fHistEt" + fHistogramNameSuffix;
 280     fHistEt = new TH1F(histname.Data(), "Total E_{T} Distribution", nbinsEt, minEt, maxEt);
 281     fHistEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
 282     fHistEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
 283
 284     histname = "fHistNeutralMult" + fHistogramNameSuffix;
 285     fHistNeutralMult = new TH1F(histname.Data(), "Neutral Multiplicity", nbinsMult, minMult, maxMult);
 286     fHistNeutralMult->GetXaxis()->SetTitle("N");
 287     fHistNeutralMult->GetYaxis()->SetTitle("Multiplicity");
 288
 289     histname = "fHistPhivsPtPos" + fHistogramNameSuffix;
 290     fHistPhivsPtPos = new TH2F(histname.Data(), "Phi vs pT of positively charged tracks hitting the calorimeter",       200, 0, 2*TMath::Pi(), nbinsPt, minPt, maxPt);
 291
 292     histname = "fHistPhivsPtNeg" + fHistogramNameSuffix;
 293     fHistPhivsPtNeg = new TH2F(histname.Data(), "Phi vs pT of negatively charged tracks hitting the calorimeter",       200, 0, 2*TMath::Pi(), nbinsPt, minPt, maxPt);
 294
 295     histname = "fCutFlow" + fHistogramNameSuffix;
 296     fCutFlow = new TH1I(histname, histname, 20, -0.5, 19.5);
 297 }
 298
 299 TH2F* AliAnalysisEt::CreateEtaEHisto2D(TString name, TString title, TString ztitle)
 300 {     //creates a 2-d histogram in eta and E and adds it to the list of histograms to be saved
 301         TString histoname   = name + fHistogramNameSuffix;
 302
 303         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfEBins, fgEAxis, fgnumOfEtaBins, fgEtaAxis);
 304         histo->SetYTitle("#eta");
 305         histo->SetXTitle("E (GeV)");
 306         histo->SetZTitle(ztitle.Data());
 307         histo->Sumw2();
 308
 309         return histo;
 310 }
 311
 312
 313 TH2F* AliAnalysisEt::CreateEtaPtHisto2D(TString name, TString title, TString ztitle)
 314 {     //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
 315         TString histoname   = name + fHistogramNameSuffix;
 316
 317         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfPtBins, fgPtAxis, fgnumOfEtaBins, fgEtaAxis);
 318         histo->SetYTitle("#eta");
 319         histo->SetXTitle("p_{T}");
 320         histo->SetZTitle(ztitle.Data());
 321         histo->Sumw2();
 322
 323         return histo;
 324 }
 325
 326 TH2F* AliAnalysisEt::CreateEtaEtHisto2D(TString name, TString title, TString ztitle)
 327 {     //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
 328         TString histoname   = name + fHistogramNameSuffix;
 329
 330         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfEBins, fgEAxis, fgnumOfEtaBins, fgEtaAxis);
 331         histo->SetYTitle("#eta");
 332         histo->SetXTitle("E_{T}");
 333         histo->SetZTitle(ztitle.Data());
 334         histo->Sumw2();
 335
 336         return histo;
 337 }
 338
 339 TH2F* AliAnalysisEt::CreateResEHisto2D(TString name, TString title, TString ztitle)
 340 {     //creates a 2-d histogram in eta and E and adds it to the list of histograms to be saved
 341         TString histoname   = name + fHistogramNameSuffix;
 342
 343         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfEBins, fgEAxis, fgNumOfRBins, fgRAxis);
 344         histo->SetYTitle("R");
 345         histo->SetXTitle("E (GeV)");
 346         histo->SetZTitle(ztitle.Data());
 347         histo->Sumw2();
 348
 349         return histo;
 350 }
 351
 352
 353 TH2F* AliAnalysisEt::CreateResPtHisto2D(TString name, TString title, TString ztitle)
 354 {     //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
 355         TString histoname   = name + fHistogramNameSuffix;
 356
 357         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfPtBins, fgPtAxis, fgNumOfRBins, fgRAxis);
 358         histo->SetYTitle("R");
 359         histo->SetXTitle("p_{T}");
 360         histo->SetZTitle(ztitle.Data());
 361         histo->Sumw2();
 362
 363         return histo;
 364 }
 365
 366 THnSparseF* AliAnalysisEt::CreateClusterHistoSparse(TString name, TString title)
 367 {     //creates a 2D sparse histogram
 368         TString histoname   = name + fHistogramNameSuffix;
 369
 370         Int_t nBins[4] = {200,200,240,20};
 371         Double_t min[4] = {0.,-1.,70.,0.5};
 372         Double_t max[4] = {50.,1.,190,20.5};
 373
 374         THnSparseF *histo = new THnSparseF(histoname.Data(),title.Data(),4,nBins,min, max);
 375     histo->GetAxis(0)->SetTitle("E");
 376     histo->GetAxis(1)->SetTitle("#eta_{cluster}");
 377     histo->GetAxis(2)->SetTitle("#phi_{cluster}");
 378     histo->GetAxis(3)->SetTitle("n_{cells}");
 379         histo->Sumw2();
 380
 381         return histo;
 382 }
 383
 384 THnSparseF* AliAnalysisEt::CreateNeutralPartHistoSparse(TString name, TString title)
 385 {     //creates a sparse neutral particle histogram
 386         TString histoname   = name + fHistogramNameSuffix;
 387
 388         Int_t nBins[6] = {20,200,200,200,240,20};
 389         Double_t min[6] = {-1.,0.,0.,-1.,70.,0.5};
 390         Double_t max[6] = {1.,50.,50.,1.,190,20.5};
 391
 392         THnSparseF *histo = new THnSparseF(histoname.Data(),title.Data(),6,nBins,min, max);
 393     histo->GetAxis(0)->SetTitle("#eta");
 394     histo->GetAxis(1)->SetTitle("p_{T}");
 395     histo->GetAxis(2)->SetTitle("E");
 396     histo->GetAxis(3)->SetTitle("#eta_{cluster}");
 397     histo->GetAxis(4)->SetTitle("#phi_{cluster}");
 398     histo->GetAxis(5)->SetTitle("n_{cells}");
 399         histo->Sumw2();
 400
 401         return histo;
 402 }
 403
 404 THnSparseF* AliAnalysisEt::CreateChargedPartHistoSparse(TString name, TString title)
 405 {     //creates a sparse charged particle histogram
 406         TString histoname   = name + fHistogramNameSuffix;
 407
 408         Int_t nBins[7] = {20,200,200,200,240,20,100};
 409         Double_t min[7] = {-1.,0.,0.,-1.,70.,0.5,0.};
 410         Double_t max[7] = {1.,50.,50.,1.,190,20.5,0.1};
 411
 412         THnSparseF *histo = new THnSparseF(histoname.Data(),title.Data(),7,nBins,min, max);
 413     histo->GetAxis(0)->SetTitle("#eta");
 414     histo->GetAxis(1)->SetTitle("p_{T}");
 415     histo->GetAxis(2)->SetTitle("E");
 416     histo->GetAxis(3)->SetTitle("#eta_{cluster}");
 417     histo->GetAxis(4)->SetTitle("#phi_{cluster}");
 418     histo->GetAxis(5)->SetTitle("n_{cells}");
 419     histo->GetAxis(6)->SetTitle("R_{match}");
 420         histo->Sumw2();
 421
 422         return histo;
 423 }
 424
 425
 426 void AliAnalysisEt::CreateTrees()
 427 { // create tree..
 428   TString treename = "fEventSummaryTree" + fHistogramNameSuffix;
 429   if(fCuts->GetHistMakeTree())
 430   {
 431
 432     fEventSummaryTree = new TTree(treename, treename);
 433     fEventSummaryTree->Branch("fTotEt",&fTotEt,"fTotEt/D");
 434     fEventSummaryTree->Branch("fTotEtAcc",&fTotEtAcc,"fTotEtAcc/D");
 435     fEventSummaryTree->Branch("fTotNeutralEt",&fTotNeutralEt,"fTotNeutralEt/D");
 436     fEventSummaryTree->Branch("fTotNeutralEtAcc",&fTotNeutralEtAcc,"fTotNeutralEtAcc/D");
 437     fEventSummaryTree->Branch("fTotChargedEt",&fTotChargedEt,"fTotChargedEt/D");
 438     fEventSummaryTree->Branch("fTotChargedEtAcc",&fTotChargedEtAcc,"fTotChargedEtAcc/D");
 439     fEventSummaryTree->Branch("fMultiplicity",&fMultiplicity,"fMultiplicity/I");
 440     fEventSummaryTree->Branch("fChargedMultiplicity",&fChargedMultiplicity,"fChargedMultiplicity/I");
 441     fEventSummaryTree->Branch("fNeutralMultiplicity",&fNeutralMultiplicity,"fNeutralMultiplicity/I");
 442     fEventSummaryTree->Branch("fCentClass",&fCentClass,"fCentClass/I");
 443     fEventSummaryTree->Branch("fChargedEnergyRemoved", &fChargedEnergyRemoved, "fChargedEnergyRemoved/D");
 444     fEventSummaryTree->Branch("fNeutralEnergyRemoved", &fNeutralEnergyRemoved, "fNeutralEnergyRemoved/D");
 445     fEventSummaryTree->Branch("fGammaEnergyAdded", &fGammaEnergyAdded, "fGammaEnergyAdded/D");
 446
 447
 448     fEventSummaryTree->Branch("fProtonEt",&fProtonEt,"fProtonEt/D");
 449     fEventSummaryTree->Branch("fAntiProtonEt",&fAntiProtonEt,"fAntiProtonEt/D");
 450
 451     fEventSummaryTree->Branch("fNeutronEt",&fNeutronEt,"fNeutronEt/D");
 452     fEventSummaryTree->Branch("fAntiNeutronEt",&fAntiNeutronEt,"fAntiNeutronEt/D");
 453
 454     fEventSummaryTree->Branch("fPi0Et",&fPi0Et,"fPi0Et/D");
 455     fEventSummaryTree->Branch("fPiPlusEt",&fPiPlusEt,"fPiPlusEt/D");
 456     fEventSummaryTree->Branch("fPiMinusEt",&fPiMinusEt,"fPiMinusEt/D");
 457
 458     fEventSummaryTree->Branch("fKPlusEt",&fKPlusEt,"fKPlusEt/D");
 459     fEventSummaryTree->Branch("fKMinusEt",&fKMinusEt,"fKMinusEt/D");
 460     fEventSummaryTree->Branch("fK0sEt",&fK0sEt,"fK0sEt/D");
 461     fEventSummaryTree->Branch("fK0lEt",&fK0lEt,"fK0lEt/D");
 462
 463     fEventSummaryTree->Branch("fMuMinusEt",&fMuMinusEt,"fMuMinusEt/D");
 464     fEventSummaryTree->Branch("fMuPlusEt",&fMuPlusEt,"fMuPlusEt/D");
 465
 466     fEventSummaryTree->Branch("fEMinusEt",&fEMinusEt,"fEMinusEt/D");
 467     fEventSummaryTree->Branch("fEPlusEt",&fEPlusEt,"fEPlusEt/D");
 468
 469     fEventSummaryTree->Branch("fGammaEt", &fGammaEt, "fGammaEt/D");
 470
 471     fEventSummaryTree->Branch("fProtonRemovedEt",&fProtonRemovedEt,"fProtonRemovedEt/D");
 472     fEventSummaryTree->Branch("fAntiProtonRemovedEt",&fAntiProtonRemovedEt,"fAntiProtonRemovedEt/D");
 473
 474     fEventSummaryTree->Branch("fNeutronRemovedEt",&fNeutronRemovedEt,"fNeutronRemovedEt/D");
 475     fEventSummaryTree->Branch("fAntiNeutronRemovedEt",&fAntiNeutronRemovedEt,"fAntiNeutronRemovedEt/D");
 476
 477     fEventSummaryTree->Branch("fPi0RemovedEt",&fPi0RemovedEt,"fPi0RemovedEt/D");
 478     fEventSummaryTree->Branch("fPiPlusRemovedEt",&fPiPlusRemovedEt,"fPiPlusRemovedEt/D");
 479     fEventSummaryTree->Branch("fPiMinusRemovedEt",&fPiMinusRemovedEt,"fPiMinusRemovedEt/D");
 480
 481     fEventSummaryTree->Branch("fKPlusRemovedEt",&fKPlusRemovedEt,"fKPlusRemovedEt/D");
 482     fEventSummaryTree->Branch("fKMinusRemovedEt",&fKMinusEt,"fKMinusRemovedEt/D");
 483     fEventSummaryTree->Branch("fK0sRemovedEt",&fK0sEt,"fK0sRemovedEt/D");
 484     fEventSummaryTree->Branch("fK0lRemovedEt",&fK0lRemovedEt,"fK0lRemovedEt/D");
 485
 486     fEventSummaryTree->Branch("fMuMinusRemovedEt",&fMuMinusRemovedEt,"fMuMinusRemovedEt/D");
 487     fEventSummaryTree->Branch("fMuPlusRemovedEt",&fMuPlusRemovedEt,"fMuPlusRemovedEt/D");
 488
 489     fEventSummaryTree->Branch("fEMinusRemovedEt",&fEMinusRemovedEt,"fEMinusRemovedEt/D");
 490     fEventSummaryTree->Branch("fEPlusRemovedEt",&fEPlusRemovedEt,"fEPlusRemovedEt/D");
 491
 492     fEventSummaryTree->Branch("fGammaRemovedEt", &fGammaRemovedEt, "fGammaEtRemoved/D");
 493
 494     fEventSummaryTree->Branch("fProtonMult",&fProtonMult,"fProtonMult/D");
 495     fEventSummaryTree->Branch("fAntiProtonMult",&fAntiProtonMult,"fAntiProtonMult/D");
 496
 497     fEventSummaryTree->Branch("fNeutronMult",&fNeutronMult,"fNeutronMult/D");
 498     fEventSummaryTree->Branch("fAntiNeutronMult",&fAntiNeutronMult,"fAntiNeutronMult/D");
 499
 500     fEventSummaryTree->Branch("fPi0Mult",&fPi0Mult,"fPi0Mult/D");
 501     fEventSummaryTree->Branch("fPiPlusMult",&fPiPlusMult,"fPiPlusMult/D");
 502     fEventSummaryTree->Branch("fPiMinusMult",&fPiMinusMult,"fPiMinusMult/D");
 503
 504     fEventSummaryTree->Branch("fKPlusMult",&fKPlusMult,"fKPlusMult/D");
 505     fEventSummaryTree->Branch("fKMinusMult",&fKMinusMult,"fKMinusMult/D");
 506     fEventSummaryTree->Branch("fK0sMult",&fK0sMult,"fK0sMult/D");
 507     fEventSummaryTree->Branch("fK0lMult",&fK0lMult,"fK0lMult/D");
 508
 509     fEventSummaryTree->Branch("fMuMinusMult",&fMuMinusMult,"fMuMinusMult/D");
 510     fEventSummaryTree->Branch("fMuPlusMult",&fMuPlusMult,"fMuPlusMult/D");
 511
 512     fEventSummaryTree->Branch("fEMinusMult",&fEMinusMult,"fEMinusMult/D");
 513     fEventSummaryTree->Branch("fEPlusMult",&fEPlusMult,"fEPlusMult/D");
 514
 515     fEventSummaryTree->Branch("fGammaMult", &fGammaMult, "fGammaMult/D");
 516
 517     fEventSummaryTree->Branch("fProtonRemovedMult",&fProtonRemovedMult,"fProtonRemovedMult/D");
 518     fEventSummaryTree->Branch("fAntiProtonRemovedMult",&fAntiProtonRemovedMult,"fAntiProtonRemovedMult/D");
 519
 520     fEventSummaryTree->Branch("fNeutronRemovedMult",&fNeutronRemovedMult,"fNeutronRemovedMult/D");
 521     fEventSummaryTree->Branch("fAntiNeutronRemovedMult",&fAntiNeutronRemovedMult,"fAntiNeutronRemovedMult/D");
 522
 523     fEventSummaryTree->Branch("fPi0RemovedMult",&fPi0RemovedMult,"fPi0RemovedMult/D");
 524     fEventSummaryTree->Branch("fPiPlusRemovedMult",&fPiPlusRemovedMult,"fPiPlusRemovedMult/D");
 525     fEventSummaryTree->Branch("fPiMinusRemovedMult",&fPiMinusRemovedMult,"fPiMinusRemovedMult/D");
 526
 527     fEventSummaryTree->Branch("fKPlusRemovedMult",&fKPlusRemovedMult,"fKPlusRemovedMult/D");
 528     fEventSummaryTree->Branch("fKMinusRemovedMult",&fKMinusMult,"fKMinusRemovedMult/D");
 529     fEventSummaryTree->Branch("fK0sRemovedMult",&fK0sMult,"fK0sRemovedMult/D");
 530     fEventSummaryTree->Branch("fK0lRemovedMult",&fK0lRemovedMult,"fK0lRemovedMult/D");
 531
 532     fEventSummaryTree->Branch("fMuMinusRemovedMult",&fMuMinusRemovedMult,"fMuMinusRemovedMult/D");
 533     fEventSummaryTree->Branch("fMuPlusRemovedMult",&fMuPlusRemovedMult,"fMuPlusRemovedMult/D");
 534
 535     fEventSummaryTree->Branch("fEMinusRemovedMult",&fEMinusRemovedMult,"fEMinusRemovedMult/D");
 536     fEventSummaryTree->Branch("fEPlusRemovedMult",&fEPlusRemovedMult,"fEPlusRemovedMult/D");
 537
 538     fEventSummaryTree->Branch("fGammaRemovedMult", &fGammaRemovedMult, "fGammaMultRemoved/D");
 539
 540
 541
 542   }
 543
 544   if(fCuts->GetHistMakeTreeDeposit())
 545   {
 546     treename = "fTreeDeposit" + fHistogramNameSuffix;
 547     fDepositTree = new TTree(treename, treename);
 548
 549     fDepositTree->Branch("fEnergyDeposited", &fEnergyDeposited, "fEnergyDeposited/F");
 550     fDepositTree->Branch("fMomentumTPC", &fMomentumTPC, "fMomentumTPC/F");
 551     fDepositTree->Branch("fCharge", &fCharge, "fCharge/S");
 552     fDepositTree->Branch("fParticlePid", &fParticlePid, "fParticlePid/S");
 553     fDepositTree->Branch("fPidProb", &fPidProb, "fPidProb/F");
 554     fDepositTree->Branch("fTrackPassedCut", &fTrackPassedCut, "fTrackPassedCut/B");
 555   }
 556
 557   return;
 558 }
 559 void AliAnalysisEt::FillHistograms()
 560 { // fill histograms..
 561     fHistEt->Fill(fTotEt);
 562
 563     fHistNeutralMult->Fill(fNeutralMultiplicity);
 564
 565     if (fCuts) {
 566       if (fCuts->GetHistMakeTree()) {
 567         fEventSummaryTree->Fill();
 568       }
 569     }
 570     if(fCuts)
 571     {
 572       if(fCuts->GetHistMakeTreeDeposit())
 573       {
 574         fDepositTree->Fill();
 575       }
 576     }
 577
 578 }
 579
 580 Int_t AliAnalysisEt::AnalyseEvent(AliVEvent *event)
 581 { //this line is basically here to eliminate a compiler warning that event is not used.  Making it a virtual function did not work with the plugin.
 582   AliAnalysisEtCommon::AnalyseEvent(event);
 583   //fSelector->SetEvent(event);
 584   ResetEventValues();
 585   return 0;
 586 }
 587
 588 void AliAnalysisEt::ResetEventValues()
 589 { // clear
 590   AliAnalysisEtCommon::ResetEventValues();
 591   fTotEt = 0;
 592   fTotEtAcc = 0;
 593   fTotNeutralEt = 0;
 594   fTotChargedEt = 0;
 595   fTotNeutralEtAcc = 0;
 596   fTotChargedEtAcc  = 0;
 597   fMultiplicity = 0;
 598   fChargedMultiplicity = 0;
 599   fNeutralMultiplicity = 0;
 600
 601   fProtonEt = 0;
 602   fAntiProtonEt = 0;
 603
 604   fNeutronEt = 0;
 605   fAntiNeutronEt = 0;
 606
 607   fPi0Et = 0;
 608   fPiPlusEt = 0;
 609   fPiMinusEt = 0;
 610
 611   fKPlusEt = 0;
 612   fKMinusEt = 0;
 613   fK0sEt = 0;
 614   fK0lEt = 0;
 615
 616   fMuMinusEt = 0;
 617   fMuPlusEt = 0;
 618
 619   fEMinusEt = 0;
 620   fEPlusEt = 0;
 621
 622   fGammaEt = 0;
 623
 624   fProtonRemovedEt = 0;
 625   fAntiProtonRemovedEt = 0;
 626
 627   fNeutronRemovedEt = 0;
 628   fAntiNeutronRemovedEt = 0;
 629
 630   fPi0RemovedEt = 0;
 631   fPiPlusRemovedEt = 0;
 632   fPiMinusRemovedEt = 0;
 633
 634   fKPlusRemovedEt = 0;
 635   fKMinusRemovedEt = 0;
 636   fK0sRemovedEt = 0;
 637   fK0lRemovedEt = 0;
 638
 639   fMuMinusRemovedEt = 0;
 640   fMuPlusRemovedEt = 0;
 641
 642   fEMinusRemovedEt = 0;
 643   fEPlusRemovedEt = 0;
 644
 645   fGammaRemovedEt = 0;
 646
 647
 648   fProtonMult = 0;
 649   fAntiProtonMult = 0;
 650
 651   fNeutronMult = 0;
 652   fAntiNeutronMult = 0;
 653
 654   fPi0Mult = 0;
 655   fPiPlusMult = 0;
 656   fPiMinusMult = 0;
 657
 658   fKPlusMult = 0;
 659   fKMinusMult = 0;
 660   fK0sMult = 0;
 661   fK0lMult = 0;
 662
 663   fMuMinusMult = 0;
 664   fMuPlusMult = 0;
 665
 666   fEMinusMult = 0;
 667   fEPlusMult = 0;
 668
 669   fGammaMult = 0;
 670
 671   fProtonRemovedMult = 0;
 672   fAntiProtonRemovedMult = 0;
 673
 674   fNeutronRemovedMult = 0;
 675   fAntiNeutronRemovedMult = 0;
 676
 677   fPi0RemovedMult = 0;
 678   fPiPlusRemovedMult = 0;
 679   fPiMinusRemovedMult = 0;
 680
 681   fKPlusRemovedMult = 0;
 682   fKMinusRemovedMult = 0;
 683   fK0sRemovedMult = 0;
 684   fK0lRemovedMult = 0;
 685
 686   fMuMinusRemovedMult = 0;
 687   fMuPlusRemovedMult = 0;
 688
 689   fEMinusRemovedMult = 0;
 690   fEPlusRemovedMult = 0;
 691
 692   fGammaRemovedMult = 0;
 693
 694   return;
 695 }
 696
 697 Int_t AliAnalysisEt::ReadCorrections(TString filename)
 698 {
 699   TFile *f = TFile::Open(filename, "READ");
 700   if(f)
 701   {
 702     TString det = "Phos";
 703     if(fHistogramNameSuffix.Contains("Emcal"))
 704     {
 705       det = "Emcal";
 706     }
 707     cout<<"Histo name suffix "<<fHistogramNameSuffix<<endl;
 708     TString name = "TmCorrections" + det;
 709     std::cout << name << std::endl;
 710     fTmCorrections = dynamic_cast<AliAnalysisEtTrackMatchCorrections*>(f->Get(name));
 711     if(!fTmCorrections)
 712     {
 713       cout<<"No corrections with name "<<name<<endl;
 714       Printf("Could not load TM corrections");
 715       return -1;
 716     }
 717     name = "ReCorrections" + det;
 718     fReCorrections = dynamic_cast<AliAnalysisEtRecEffCorrection*>(f->Get(name));
 719     if(!fReCorrections)
 720     {
 721       Printf("Could not load rec eff corrections");
 722       return -1;
 723     }
 724     return 0;
 725   }
 726  return -1;
 727 }
 728
 729 Double_t AliAnalysisEt::CorrectForReconstructionEfficiency(const AliESDCaloCluster& cluster, Int_t cent)
 730 {
 731   Float_t pos[3];
 732   cluster.GetPosition(pos);
 733   TVector3 cp(pos);
 734   Double_t corrEnergy = fReCorrections->CorrectedEnergy(cluster.E(), cent);
 735
 736   //std::cout << "Original energy: " << cluster.E() << ", corrected energy: " << corrEnergy << std::endl;
 737   return TMath::Sin(cp.Theta())*corrEnergy;
 738 }
 739
 740
 741 Double_t AliAnalysisEt::CorrectForReconstructionEfficiency(const AliESDCaloCluster& cluster, Float_t eReco, Int_t cent)
 742 {
 743   Float_t pos[3];
 744   cluster.GetPosition(pos);
 745   TVector3 cp(pos);
 746   Double_t corrEnergy = fReCorrections->CorrectedEnergy(eReco, cent);
 747
 748   //std::cout << "Original energy: " << cluster.E() << ", corrected energy: " << corrEnergy << std::endl;
 749   return TMath::Sin(cp.Theta())*corrEnergy;
 750 }
 751