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
 166 {}
 167
 168 AliAnalysisEt::~AliAnalysisEt()
 169 {//Destructor
 170   delete fTmCorrections;
 171   delete fReCorrections;
 172   if(fDepositTree){
 173     fDepositTree->Clear();
 174     delete fDepositTree; // optional TTree
 175   }
 176   if(fEventSummaryTree)
 177   {
 178     fEventSummaryTree->Clear();
 179     delete fEventSummaryTree;
 180   }
 181   if(fAcceptedTree)
 182   {
 183     fAcceptedTree->Clear();
 184     delete fAcceptedTree;
 185   }
 186   delete fHistEt; //Et spectrum
 187   delete fHistNeutralMult; //Neutral multiplicity
 188   delete fHistPhivsPtPos; //phi vs pT plot for positive tracks
 189   delete fHistPhivsPtNeg; //phi vs pT Moplot for negative tracks
 190   //delete fCentrality;//this code does not actually own AliCentrality so we don't have to worry about deleting it...  we just borrow it...
 191   delete fCutFlow;
 192   delete fSelector;
 193   delete fPIDResponse;
 194 }
 195
 196 void AliAnalysisEt::FillOutputList(TList *list)
 197 { // histograms to be added to output
 198     list->Add(fHistEt);
 199     list->Add(fHistNeutralMult);
 200
 201     list->Add(fHistPhivsPtPos);
 202     list->Add(fHistPhivsPtNeg);
 203
 204     if (fCuts) {
 205       if (fCuts->GetHistMakeTree()) {
 206         //list->Add(fTree);
 207         list->Add(fEventSummaryTree);
 208       }
 209       if (fCuts->GetHistMakeTreeDeposit()) {
 210         list->Add(fDepositTree);
 211       }
 212     }
 213
 214     list->Add(fCutFlow);
 215
 216     AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 217     if (!man) {
 218       AliFatal("Analysis manager needed");
 219       return;
 220     }
 221
 222     AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(man->GetInputEventHandler());
 223     if (!inputHandler) {
 224       AliFatal("Input handler needed");
 225       return;
 226     }
 227
 228     //pid response object
 229     fPIDResponse=inputHandler->GetPIDResponse();
 230     if (!fPIDResponse) AliError("PIDResponse object was not created");
 231
 232 }
 233
 234 void AliAnalysisEt::Init()
 235 {// clear variables, set up cuts and PDG info
 236   AliAnalysisEtCommon::Init();
 237   if(ReadCorrections("calocorrections.root") != 0)
 238   {
 239     // Shouldn't do this, why oh why are exceptions not allowed?
 240     exit(-1);
 241   }
 242   ResetEventValues();
 243
 244 }
 245
 246 void AliAnalysisEt::CreateHistograms()
 247 { // create histograms..
 248   // histogram binning for E_T, p_T and Multiplicity: defaults for p+p
 249   Int_t nbinsEt = 10000;
 250   Double_t minEt = 0.0;
 251   Double_t maxEt = 1000;
 252   Int_t nbinsPt = 200;
 253   Double_t minPt = 0;
 254   Double_t maxPt = 20;
 255   Int_t nbinsMult = 200;
 256   Double_t minMult = -0.5; // offset -0.5 to have integer bins centered around 0
 257   Double_t maxMult = nbinsMult + minMult; // 1 bin per integer value
 258
 259   // see if we should change histogram limits etc, and possibly create a tree
 260   if (fCuts) {
 261     if (fCuts->GetHistMakeTree()) {
 262       CreateTrees();
 263     }
 264
 265     nbinsMult = fCuts->GetHistNbinsMult();
 266     minMult = fCuts->GetHistMinMult();
 267     maxMult = fCuts->GetHistMaxMult();
 268
 269     nbinsEt = fCuts->GetHistNbinsTotEt();
 270     minEt = fCuts->GetHistMinTotEt();
 271     maxEt = fCuts->GetHistMaxTotEt();
 272
 273     nbinsPt = fCuts->GetHistNbinsParticlePt();
 274     minPt = fCuts->GetHistMinParticlePt();
 275     maxPt = fCuts->GetHistMaxParticlePt();
 276   }
 277
 278     TString histname = "fHistEt" + fHistogramNameSuffix;
 279     fHistEt = new TH1F(histname.Data(), "Total E_{T} Distribution", nbinsEt, minEt, maxEt);
 280     fHistEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
 281     fHistEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
 282
 283     histname = "fHistNeutralMult" + fHistogramNameSuffix;
 284     fHistNeutralMult = new TH1F(histname.Data(), "Neutral Multiplicity", nbinsMult, minMult, maxMult);
 285     fHistNeutralMult->GetXaxis()->SetTitle("N");
 286     fHistNeutralMult->GetYaxis()->SetTitle("Multiplicity");
 287
 288     histname = "fHistPhivsPtPos" + fHistogramNameSuffix;
 289     fHistPhivsPtPos = new TH2F(histname.Data(), "Phi vs pT of positively charged tracks hitting the calorimeter",       200, 0, 2*TMath::Pi(), nbinsPt, minPt, maxPt);
 290
 291     histname = "fHistPhivsPtNeg" + fHistogramNameSuffix;
 292     fHistPhivsPtNeg = new TH2F(histname.Data(), "Phi vs pT of negatively charged tracks hitting the calorimeter",       200, 0, 2*TMath::Pi(), nbinsPt, minPt, maxPt);
 293
 294     histname = "fCutFlow" + fHistogramNameSuffix;
 295     fCutFlow = new TH1I(histname, histname, 20, -0.5, 19.5);
 296 }
 297
 298 TH2F* AliAnalysisEt::CreateEtaEHisto2D(TString name, TString title, TString ztitle)
 299 {     //creates a 2-d histogram in eta and E and adds it to the list of histograms to be saved
 300         TString histoname   = name + fHistogramNameSuffix;
 301
 302         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfEBins, fgEAxis, fgnumOfEtaBins, fgEtaAxis);
 303         histo->SetYTitle("#eta");
 304         histo->SetXTitle("E (GeV)");
 305         histo->SetZTitle(ztitle.Data());
 306         histo->Sumw2();
 307
 308         return histo;
 309 }
 310
 311
 312 TH2F* AliAnalysisEt::CreateEtaPtHisto2D(TString name, TString title, TString ztitle)
 313 {     //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
 314         TString histoname   = name + fHistogramNameSuffix;
 315
 316         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfPtBins, fgPtAxis, fgnumOfEtaBins, fgEtaAxis);
 317         histo->SetYTitle("#eta");
 318         histo->SetXTitle("p_{T}");
 319         histo->SetZTitle(ztitle.Data());
 320         histo->Sumw2();
 321
 322         return histo;
 323 }
 324
 325 TH2F* AliAnalysisEt::CreateEtaEtHisto2D(TString name, TString title, TString ztitle)
 326 {     //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
 327         TString histoname   = name + fHistogramNameSuffix;
 328
 329         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfEBins, fgEAxis, fgnumOfEtaBins, fgEtaAxis);
 330         histo->SetYTitle("#eta");
 331         histo->SetXTitle("E_{T}");
 332         histo->SetZTitle(ztitle.Data());
 333         histo->Sumw2();
 334
 335         return histo;
 336 }
 337
 338 TH2F* AliAnalysisEt::CreateResEHisto2D(TString name, TString title, TString ztitle)
 339 {     //creates a 2-d histogram in eta and E and adds it to the list of histograms to be saved
 340         TString histoname   = name + fHistogramNameSuffix;
 341
 342         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfEBins, fgEAxis, fgNumOfRBins, fgRAxis);
 343         histo->SetYTitle("R");
 344         histo->SetXTitle("E (GeV)");
 345         histo->SetZTitle(ztitle.Data());
 346         histo->Sumw2();
 347
 348         return histo;
 349 }
 350
 351
 352 TH2F* AliAnalysisEt::CreateResPtHisto2D(TString name, TString title, TString ztitle)
 353 {     //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
 354         TString histoname   = name + fHistogramNameSuffix;
 355
 356         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfPtBins, fgPtAxis, fgNumOfRBins, fgRAxis);
 357         histo->SetYTitle("R");
 358         histo->SetXTitle("p_{T}");
 359         histo->SetZTitle(ztitle.Data());
 360         histo->Sumw2();
 361
 362         return histo;
 363 }
 364
 365 THnSparseF* AliAnalysisEt::CreateClusterHistoSparse(TString name, TString title)
 366 {     //creates a 2D sparse histogram
 367         TString histoname   = name + fHistogramNameSuffix;
 368
 369         Int_t nBins[4] = {200,200,240,20};
 370         Double_t min[4] = {0.,-1.,70.,0.5};
 371         Double_t max[4] = {50.,1.,190,20.5};
 372
 373         THnSparseF *histo = new THnSparseF(histoname.Data(),title.Data(),4,nBins,min, max);
 374     histo->GetAxis(0)->SetTitle("E");
 375     histo->GetAxis(1)->SetTitle("#eta_{cluster}");
 376     histo->GetAxis(2)->SetTitle("#phi_{cluster}");
 377     histo->GetAxis(3)->SetTitle("n_{cells}");
 378         histo->Sumw2();
 379
 380         return histo;
 381 }
 382
 383 THnSparseF* AliAnalysisEt::CreateNeutralPartHistoSparse(TString name, TString title)
 384 {     //creates a sparse neutral particle histogram
 385         TString histoname   = name + fHistogramNameSuffix;
 386
 387         Int_t nBins[6] = {20,200,200,200,240,20};
 388         Double_t min[6] = {-1.,0.,0.,-1.,70.,0.5};
 389         Double_t max[6] = {1.,50.,50.,1.,190,20.5};
 390
 391         THnSparseF *histo = new THnSparseF(histoname.Data(),title.Data(),6,nBins,min, max);
 392     histo->GetAxis(0)->SetTitle("#eta");
 393     histo->GetAxis(1)->SetTitle("p_{T}");
 394     histo->GetAxis(2)->SetTitle("E");
 395     histo->GetAxis(3)->SetTitle("#eta_{cluster}");
 396     histo->GetAxis(4)->SetTitle("#phi_{cluster}");
 397     histo->GetAxis(5)->SetTitle("n_{cells}");
 398         histo->Sumw2();
 399
 400         return histo;
 401 }
 402
 403 THnSparseF* AliAnalysisEt::CreateChargedPartHistoSparse(TString name, TString title)
 404 {     //creates a sparse charged particle histogram
 405         TString histoname   = name + fHistogramNameSuffix;
 406
 407         Int_t nBins[7] = {20,200,200,200,240,20,100};
 408         Double_t min[7] = {-1.,0.,0.,-1.,70.,0.5,0.};
 409         Double_t max[7] = {1.,50.,50.,1.,190,20.5,0.1};
 410
 411         THnSparseF *histo = new THnSparseF(histoname.Data(),title.Data(),7,nBins,min, max);
 412     histo->GetAxis(0)->SetTitle("#eta");
 413     histo->GetAxis(1)->SetTitle("p_{T}");
 414     histo->GetAxis(2)->SetTitle("E");
 415     histo->GetAxis(3)->SetTitle("#eta_{cluster}");
 416     histo->GetAxis(4)->SetTitle("#phi_{cluster}");
 417     histo->GetAxis(5)->SetTitle("n_{cells}");
 418     histo->GetAxis(6)->SetTitle("R_{match}");
 419         histo->Sumw2();
 420
 421         return histo;
 422 }
 423
 424
 425 void AliAnalysisEt::CreateTrees()
 426 { // create tree..
 427   TString treename = "fEventSummaryTree" + fHistogramNameSuffix;
 428   if(fCuts->GetHistMakeTree())
 429   {
 430
 431     fEventSummaryTree = new TTree(treename, treename);
 432     fEventSummaryTree->Branch("fTotEt",&fTotEt,"fTotEt/D");
 433     fEventSummaryTree->Branch("fTotEtAcc",&fTotEtAcc,"fTotEtAcc/D");
 434     fEventSummaryTree->Branch("fTotNeutralEt",&fTotNeutralEt,"fTotNeutralEt/D");
 435     fEventSummaryTree->Branch("fTotNeutralEtAcc",&fTotNeutralEtAcc,"fTotNeutralEtAcc/D");
 436     fEventSummaryTree->Branch("fTotChargedEt",&fTotChargedEt,"fTotChargedEt/D");
 437     fEventSummaryTree->Branch("fTotChargedEtAcc",&fTotChargedEtAcc,"fTotChargedEtAcc/D");
 438     fEventSummaryTree->Branch("fMultiplicity",&fMultiplicity,"fMultiplicity/I");
 439     fEventSummaryTree->Branch("fChargedMultiplicity",&fChargedMultiplicity,"fChargedMultiplicity/I");
 440     fEventSummaryTree->Branch("fNeutralMultiplicity",&fNeutralMultiplicity,"fNeutralMultiplicity/I");
 441     fEventSummaryTree->Branch("fCentClass",&fCentClass,"fCentClass/I");
 442     fEventSummaryTree->Branch("fChargedEnergyRemoved", &fChargedEnergyRemoved, "fChargedEnergyRemoved/D");
 443     fEventSummaryTree->Branch("fNeutralEnergyRemoved", &fNeutralEnergyRemoved, "fNeutralEnergyRemoved/D");
 444     fEventSummaryTree->Branch("fGammaEnergyAdded", &fGammaEnergyAdded, "fGammaEnergyAdded/D");
 445
 446
 447     fEventSummaryTree->Branch("fProtonEt",&fProtonEt,"fProtonEt/D");
 448     fEventSummaryTree->Branch("fAntiProtonEt",&fAntiProtonEt,"fAntiProtonEt/D");
 449
 450     fEventSummaryTree->Branch("fNeutronEt",&fNeutronEt,"fNeutronEt/D");
 451     fEventSummaryTree->Branch("fAntiNeutronEt",&fAntiNeutronEt,"fAntiNeutronEt/D");
 452
 453     fEventSummaryTree->Branch("fPi0Et",&fPi0Et,"fPi0Et/D");
 454     fEventSummaryTree->Branch("fPiPlusEt",&fPiPlusEt,"fPiPlusEt/D");
 455     fEventSummaryTree->Branch("fPiMinusEt",&fPiMinusEt,"fPiMinusEt/D");
 456
 457     fEventSummaryTree->Branch("fKPlusEt",&fKPlusEt,"fKPlusEt/D");
 458     fEventSummaryTree->Branch("fKMinusEt",&fKMinusEt,"fKMinusEt/D");
 459     fEventSummaryTree->Branch("fK0sEt",&fK0sEt,"fK0sEt/D");
 460     fEventSummaryTree->Branch("fK0lEt",&fK0lEt,"fK0lEt/D");
 461
 462     fEventSummaryTree->Branch("fMuMinusEt",&fMuMinusEt,"fMuMinusEt/D");
 463     fEventSummaryTree->Branch("fMuPlusEt",&fMuPlusEt,"fMuPlusEt/D");
 464
 465     fEventSummaryTree->Branch("fEMinusEt",&fEMinusEt,"fEMinusEt/D");
 466     fEventSummaryTree->Branch("fEPlusEt",&fEPlusEt,"fEPlusEt/D");
 467
 468     fEventSummaryTree->Branch("fGammaEt", &fGammaEt, "fGammaEt/D");
 469
 470     fEventSummaryTree->Branch("fProtonRemovedEt",&fProtonRemovedEt,"fProtonRemovedEt/D");
 471     fEventSummaryTree->Branch("fAntiProtonRemovedEt",&fAntiProtonRemovedEt,"fAntiProtonRemovedEt/D");
 472
 473     fEventSummaryTree->Branch("fNeutronRemovedEt",&fNeutronRemovedEt,"fNeutronRemovedEt/D");
 474     fEventSummaryTree->Branch("fAntiNeutronRemovedEt",&fAntiNeutronRemovedEt,"fAntiNeutronRemovedEt/D");
 475
 476     fEventSummaryTree->Branch("fPi0RemovedEt",&fPi0RemovedEt,"fPi0RemovedEt/D");
 477     fEventSummaryTree->Branch("fPiPlusRemovedEt",&fPiPlusRemovedEt,"fPiPlusRemovedEt/D");
 478     fEventSummaryTree->Branch("fPiMinusRemovedEt",&fPiMinusRemovedEt,"fPiMinusRemovedEt/D");
 479
 480     fEventSummaryTree->Branch("fKPlusRemovedEt",&fKPlusRemovedEt,"fKPlusRemovedEt/D");
 481     fEventSummaryTree->Branch("fKMinusRemovedEt",&fKMinusEt,"fKMinusRemovedEt/D");
 482     fEventSummaryTree->Branch("fK0sRemovedEt",&fK0sEt,"fK0sRemovedEt/D");
 483     fEventSummaryTree->Branch("fK0lRemovedEt",&fK0lRemovedEt,"fK0lRemovedEt/D");
 484
 485     fEventSummaryTree->Branch("fMuMinusRemovedEt",&fMuMinusRemovedEt,"fMuMinusRemovedEt/D");
 486     fEventSummaryTree->Branch("fMuPlusRemovedEt",&fMuPlusRemovedEt,"fMuPlusRemovedEt/D");
 487
 488     fEventSummaryTree->Branch("fEMinusRemovedEt",&fEMinusRemovedEt,"fEMinusRemovedEt/D");
 489     fEventSummaryTree->Branch("fEPlusRemovedEt",&fEPlusRemovedEt,"fEPlusRemovedEt/D");
 490
 491     fEventSummaryTree->Branch("fGammaRemovedEt", &fGammaRemovedEt, "fGammaEtRemoved/D");
 492
 493     fEventSummaryTree->Branch("fProtonMult",&fProtonMult,"fProtonMult/D");
 494     fEventSummaryTree->Branch("fAntiProtonMult",&fAntiProtonMult,"fAntiProtonMult/D");
 495
 496     fEventSummaryTree->Branch("fNeutronMult",&fNeutronMult,"fNeutronMult/D");
 497     fEventSummaryTree->Branch("fAntiNeutronMult",&fAntiNeutronMult,"fAntiNeutronMult/D");
 498
 499     fEventSummaryTree->Branch("fPi0Mult",&fPi0Mult,"fPi0Mult/D");
 500     fEventSummaryTree->Branch("fPiPlusMult",&fPiPlusMult,"fPiPlusMult/D");
 501     fEventSummaryTree->Branch("fPiMinusMult",&fPiMinusMult,"fPiMinusMult/D");
 502
 503     fEventSummaryTree->Branch("fKPlusMult",&fKPlusMult,"fKPlusMult/D");
 504     fEventSummaryTree->Branch("fKMinusMult",&fKMinusMult,"fKMinusMult/D");
 505     fEventSummaryTree->Branch("fK0sMult",&fK0sMult,"fK0sMult/D");
 506     fEventSummaryTree->Branch("fK0lMult",&fK0lMult,"fK0lMult/D");
 507
 508     fEventSummaryTree->Branch("fMuMinusMult",&fMuMinusMult,"fMuMinusMult/D");
 509     fEventSummaryTree->Branch("fMuPlusMult",&fMuPlusMult,"fMuPlusMult/D");
 510
 511     fEventSummaryTree->Branch("fEMinusMult",&fEMinusMult,"fEMinusMult/D");
 512     fEventSummaryTree->Branch("fEPlusMult",&fEPlusMult,"fEPlusMult/D");
 513
 514     fEventSummaryTree->Branch("fGammaMult", &fGammaMult, "fGammaMult/D");
 515
 516     fEventSummaryTree->Branch("fProtonRemovedMult",&fProtonRemovedMult,"fProtonRemovedMult/D");
 517     fEventSummaryTree->Branch("fAntiProtonRemovedMult",&fAntiProtonRemovedMult,"fAntiProtonRemovedMult/D");
 518
 519     fEventSummaryTree->Branch("fNeutronRemovedMult",&fNeutronRemovedMult,"fNeutronRemovedMult/D");
 520     fEventSummaryTree->Branch("fAntiNeutronRemovedMult",&fAntiNeutronRemovedMult,"fAntiNeutronRemovedMult/D");
 521
 522     fEventSummaryTree->Branch("fPi0RemovedMult",&fPi0RemovedMult,"fPi0RemovedMult/D");
 523     fEventSummaryTree->Branch("fPiPlusRemovedMult",&fPiPlusRemovedMult,"fPiPlusRemovedMult/D");
 524     fEventSummaryTree->Branch("fPiMinusRemovedMult",&fPiMinusRemovedMult,"fPiMinusRemovedMult/D");
 525
 526     fEventSummaryTree->Branch("fKPlusRemovedMult",&fKPlusRemovedMult,"fKPlusRemovedMult/D");
 527     fEventSummaryTree->Branch("fKMinusRemovedMult",&fKMinusMult,"fKMinusRemovedMult/D");
 528     fEventSummaryTree->Branch("fK0sRemovedMult",&fK0sMult,"fK0sRemovedMult/D");
 529     fEventSummaryTree->Branch("fK0lRemovedMult",&fK0lRemovedMult,"fK0lRemovedMult/D");
 530
 531     fEventSummaryTree->Branch("fMuMinusRemovedMult",&fMuMinusRemovedMult,"fMuMinusRemovedMult/D");
 532     fEventSummaryTree->Branch("fMuPlusRemovedMult",&fMuPlusRemovedMult,"fMuPlusRemovedMult/D");
 533
 534     fEventSummaryTree->Branch("fEMinusRemovedMult",&fEMinusRemovedMult,"fEMinusRemovedMult/D");
 535     fEventSummaryTree->Branch("fEPlusRemovedMult",&fEPlusRemovedMult,"fEPlusRemovedMult/D");
 536
 537     fEventSummaryTree->Branch("fGammaRemovedMult", &fGammaRemovedMult, "fGammaMultRemoved/D");
 538
 539
 540
 541   }
 542
 543   if(fCuts->GetHistMakeTreeDeposit())
 544   {
 545     treename = "fTreeDeposit" + fHistogramNameSuffix;
 546     fDepositTree = new TTree(treename, treename);
 547
 548     fDepositTree->Branch("fEnergyDeposited", &fEnergyDeposited, "fEnergyDeposited/F");
 549     fDepositTree->Branch("fMomentumTPC", &fMomentumTPC, "fMomentumTPC/F");
 550     fDepositTree->Branch("fCharge", &fCharge, "fCharge/S");
 551     fDepositTree->Branch("fParticlePid", &fParticlePid, "fParticlePid/S");
 552     fDepositTree->Branch("fPidProb", &fPidProb, "fPidProb/F");
 553     fDepositTree->Branch("fTrackPassedCut", &fTrackPassedCut, "fTrackPassedCut/B");
 554   }
 555
 556   return;
 557 }
 558 void AliAnalysisEt::FillHistograms()
 559 { // fill histograms..
 560     fHistEt->Fill(fTotEt);
 561
 562     fHistNeutralMult->Fill(fNeutralMultiplicity);
 563
 564     if (fCuts) {
 565       if (fCuts->GetHistMakeTree()) {
 566         fEventSummaryTree->Fill();
 567       }
 568     }
 569     if(fCuts)
 570     {
 571       if(fCuts->GetHistMakeTreeDeposit())
 572       {
 573         fDepositTree->Fill();
 574       }
 575     }
 576
 577 }
 578
 579 Int_t AliAnalysisEt::AnalyseEvent(AliVEvent *event)
 580 { //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.
 581   AliAnalysisEtCommon::AnalyseEvent(event);
 582   //fSelector->SetEvent(event);
 583   ResetEventValues();
 584   return 0;
 585 }
 586
 587 void AliAnalysisEt::ResetEventValues()
 588 { // clear
 589   AliAnalysisEtCommon::ResetEventValues();
 590   fTotEt = 0;
 591   fTotEtAcc = 0;
 592   fTotNeutralEt = 0;
 593   fTotChargedEt = 0;
 594   fTotNeutralEtAcc = 0;
 595   fTotChargedEtAcc  = 0;
 596   fMultiplicity = 0;
 597   fChargedMultiplicity = 0;
 598   fNeutralMultiplicity = 0;
 599
 600   fProtonEt = 0;
 601   fAntiProtonEt = 0;
 602
 603   fNeutronEt = 0;
 604   fAntiNeutronEt = 0;
 605
 606   fPi0Et = 0;
 607   fPiPlusEt = 0;
 608   fPiMinusEt = 0;
 609
 610   fKPlusEt = 0;
 611   fKMinusEt = 0;
 612   fK0sEt = 0;
 613   fK0lEt = 0;
 614
 615   fMuMinusEt = 0;
 616   fMuPlusEt = 0;
 617
 618   fEMinusEt = 0;
 619   fEPlusEt = 0;
 620
 621   fGammaEt = 0;
 622
 623   fProtonRemovedEt = 0;
 624   fAntiProtonRemovedEt = 0;
 625
 626   fNeutronRemovedEt = 0;
 627   fAntiNeutronRemovedEt = 0;
 628
 629   fPi0RemovedEt = 0;
 630   fPiPlusRemovedEt = 0;
 631   fPiMinusRemovedEt = 0;
 632
 633   fKPlusRemovedEt = 0;
 634   fKMinusRemovedEt = 0;
 635   fK0sRemovedEt = 0;
 636   fK0lRemovedEt = 0;
 637
 638   fMuMinusRemovedEt = 0;
 639   fMuPlusRemovedEt = 0;
 640
 641   fEMinusRemovedEt = 0;
 642   fEPlusRemovedEt = 0;
 643
 644   fGammaRemovedEt = 0;
 645
 646
 647   fProtonMult = 0;
 648   fAntiProtonMult = 0;
 649
 650   fNeutronMult = 0;
 651   fAntiNeutronMult = 0;
 652
 653   fPi0Mult = 0;
 654   fPiPlusMult = 0;
 655   fPiMinusMult = 0;
 656
 657   fKPlusMult = 0;
 658   fKMinusMult = 0;
 659   fK0sMult = 0;
 660   fK0lMult = 0;
 661
 662   fMuMinusMult = 0;
 663   fMuPlusMult = 0;
 664
 665   fEMinusMult = 0;
 666   fEPlusMult = 0;
 667
 668   fGammaMult = 0;
 669
 670   fProtonRemovedMult = 0;
 671   fAntiProtonRemovedMult = 0;
 672
 673   fNeutronRemovedMult = 0;
 674   fAntiNeutronRemovedMult = 0;
 675
 676   fPi0RemovedMult = 0;
 677   fPiPlusRemovedMult = 0;
 678   fPiMinusRemovedMult = 0;
 679
 680   fKPlusRemovedMult = 0;
 681   fKMinusRemovedMult = 0;
 682   fK0sRemovedMult = 0;
 683   fK0lRemovedMult = 0;
 684
 685   fMuMinusRemovedMult = 0;
 686   fMuPlusRemovedMult = 0;
 687
 688   fEMinusRemovedMult = 0;
 689   fEPlusRemovedMult = 0;
 690
 691   fGammaRemovedMult = 0;
 692
 693   return;
 694 }
 695
 696 Int_t AliAnalysisEt::ReadCorrections(TString filename)
 697 {
 698   TFile *f = TFile::Open(filename, "READ");
 699   if(f)
 700   {
 701     TString det = "Phos";
 702     if(fHistogramNameSuffix.Contains("Emcal"))
 703     {
 704       det = "Emcal";
 705     }
 706     cout<<"Histo name suffix "<<fHistogramNameSuffix<<endl;
 707     TString name = "TmCorrections" + det;
 708     std::cout << name << std::endl;
 709     fTmCorrections = dynamic_cast<AliAnalysisEtTrackMatchCorrections*>(f->Get(name));
 710     if(!fTmCorrections)
 711     {
 712       cout<<"No corrections with name "<<name<<endl;
 713       Printf("Could not load TM corrections");
 714       return -1;
 715     }
 716     name = "ReCorrections" + det;
 717     fReCorrections = dynamic_cast<AliAnalysisEtRecEffCorrection*>(f->Get(name));
 718     if(!fReCorrections)
 719     {
 720       Printf("Could not load rec eff corrections");
 721       return -1;
 722     }
 723     return 0;
 724   }
 725  return -1;
 726 }
 727
 728 Double_t AliAnalysisEt::CorrectForReconstructionEfficiency(const AliESDCaloCluster& cluster, Int_t cent)
 729 {
 730   Float_t pos[3];
 731   cluster.GetPosition(pos);
 732   TVector3 cp(pos);
 733   Double_t corrEnergy = fReCorrections->CorrectedEnergy(cluster.E(), cent);
 734
 735   //std::cout << "Original energy: " << cluster.E() << ", corrected energy: " << corrEnergy << std::endl;
 736   return TMath::Sin(cp.Theta())*corrEnergy;
 737 }
 738
 739
 740 Double_t AliAnalysisEt::CorrectForReconstructionEfficiency(const AliESDCaloCluster& cluster, Float_t eReco, Int_t cent)
 741 {
 742   Float_t pos[3];
 743   cluster.GetPosition(pos);
 744   TVector3 cp(pos);
 745   Double_t corrEnergy = fReCorrections->CorrectedEnergy(eReco, cent);
 746
 747   //std::cout << "Original energy: " << cluster.E() << ", corrected energy: " << corrEnergy << std::endl;
 748   return TMath::Sin(cp.Theta())*corrEnergy;
 749 }
 750