PWGLF/totEt/AliAnalysisEtReconstructed.cxx

   1 //_________________________________________________________________________
   2 //  Utility Class for transverse energy studies
   3 //  Base class for ESD analysis
   4 //  - reconstruction output
   5 //  implementation file
   6 //
   7 //*-- Authors: Oystein Djuvsland (Bergen), David Silvermyr (ORNL)
   8 //_________________________________________________________________________
   9
  10 #include "AliAnalysisEtReconstructed.h"
  11 #include "AliAnalysisEtCuts.h"
  12 #include "AliESDtrack.h"
  13 #include "AliEMCALTrack.h"
  14 #include "AliESDCaloCluster.h"
  15 #include "TVector3.h"
  16 #include "TGeoGlobalMagField.h"
  17 #include "AliMagF.h"
  18 #include "AliVEvent.h"
  19 #include "AliESDEvent.h"
  20 #include "AliESDtrackCuts.h"
  21 #include "AliVParticle.h"
  22 #include "TDatabasePDG.h"
  23 #include "TList.h"
  24 #include "AliESDpid.h"
  25 #include <iostream>
  26 #include "TH3F.h"
  27 #include "TH2F.h"
  28 #include "TH2I.h"
  29 #include "TH1I.h"
  30 #include "TFile.h"
  31 #include "AliAnalysisHadEtCorrections.h"
  32 #include "AliAnalysisEtSelector.h"
  33 #include "AliLog.h"
  34 #include "AliCentrality.h"
  35 #include "AliPHOSGeoUtils.h"
  36 #include "AliPHOSGeometry.h"
  37 #include "AliAnalysisEtRecEffCorrection.h"
  38
  39
  40 using namespace std;
  41
  42 ClassImp(AliAnalysisEtReconstructed);
  43
  44
  45 AliAnalysisEtReconstructed::AliAnalysisEtReconstructed() :
  46         AliAnalysisEt()
  47         ,fCorrections(0)
  48         ,fPidCut(0)
  49         ,fHistChargedPionEnergyDeposit(0)
  50         ,fHistProtonEnergyDeposit(0)
  51         ,fHistAntiProtonEnergyDeposit(0)
  52         ,fHistChargedKaonEnergyDeposit(0)
  53         ,fHistMuonEnergyDeposit(0)
  54         ,fHistRemovedEnergy(0)
  55         ,fGeomCorrection(1.0)
  56         ,fEMinCorrection(1.0/0.687)
  57         ,fRecEffCorrection(1.0)
  58         ,fClusterPosition(0)
  59         ,fClusterEnergy(0)
  60         ,fClusterEt(0)
  61         ,fHistChargedEnergyRemoved(0)
  62         ,fHistNeutralEnergyRemoved(0)
  63         ,fHistGammaEnergyAdded(0)
  64         ,fHistMatchedTracksEvspTvsMult(0)
  65         ,fHistMatchedTracksEvspTvsMultEffCorr(0)
  66         ,fHistFoundHadronsvsCent(0)
  67         ,fHistNotFoundHadronsvsCent(0)
  68         ,fHistFoundHadronsEtvsCent(0)
  69         ,fHistNotFoundHadronsEtvsCent(0)
  70         ,fHistNominalRawEt(0)
  71         ,fHistNominalNonLinHighEt(0)
  72         ,fHistNominalNonLinLowEt(0)
  73         ,fHistNominalEffHighEt(0)
  74         ,fHistNominalEffLowEt(0)
  75 {
  76
  77 }
  78
  79 AliAnalysisEtReconstructed::~AliAnalysisEtReconstructed()
  80 {//destructor
  81     delete fCorrections;
  82     delete fHistChargedPionEnergyDeposit; /** Energy deposited in calorimeter by charged pions */
  83     delete fHistProtonEnergyDeposit; /** Energy deposited in calorimeter by protons */
  84     delete fHistAntiProtonEnergyDeposit; /** Energy deposited in calorimeter by anti-protons */
  85     delete fHistChargedKaonEnergyDeposit; /** Energy deposited in calorimeter by charged kaons */
  86     delete fHistMuonEnergyDeposit; /** Energy deposited in calorimeter by muons */
  87
  88     delete fHistRemovedEnergy; // removed energy
  89     delete fClusterPosition;
  90     delete fClusterEnergy;
  91     delete fClusterEt;
  92     delete fHistChargedEnergyRemoved;
  93     delete fHistNeutralEnergyRemoved;
  94     delete fHistGammaEnergyAdded;
  95     delete fHistMatchedTracksEvspTvsMult;
  96     delete fHistMatchedTracksEvspTvsMultEffCorr;
  97     delete fHistFoundHadronsvsCent;
  98     delete fHistNotFoundHadronsvsCent;
  99     delete fHistFoundHadronsEtvsCent;
 100     delete fHistNotFoundHadronsEtvsCent;
 101     delete fHistNominalRawEt;
 102     delete fHistNominalNonLinHighEt;
 103     delete fHistNominalNonLinLowEt;
 104     delete fHistNominalEffHighEt;
 105     delete fHistNominalEffLowEt;
 106 }
 107
 108 Int_t AliAnalysisEtReconstructed::AnalyseEvent(AliVEvent* ev)
 109 {
 110
 111     //AliAnalysisEt::AnalyseEvent(ev);
 112     // analyse ESD event
 113     ResetEventValues();
 114     if (!ev) {
 115         AliFatal("ERROR: Event does not exist");
 116         return 0;
 117     }
 118
 119     AliESDEvent *event = dynamic_cast<AliESDEvent*>(ev);
 120     if (!event) {
 121         AliFatal("ERROR: ESD Event does not exist");
 122         return 0;
 123     }
 124     if(!fSelector){
 125         AliFatal("ERROR: fSelector does not exist");
 126         return 0;
 127     }
 128     fSelector->SetEvent(event);
 129
 130     Int_t cent = -1;
 131     fCentrality = event->GetCentrality();
 132     if (fCentrality && cent)
 133     {
 134         cent = fCentrality->GetCentralityClass5("V0M");
 135         fCentClass = fCentrality->GetCentralityClass5("V0M");
 136     }
 137
 138     TRefArray *caloClusters = fSelector->GetClusters();
 139     Float_t fClusterMult = caloClusters->GetEntries();
 140
 141     Float_t nominalRawEt = 0;
 142     Float_t nonlinHighRawEt = 0;
 143     Float_t nonlinLowRawEt = 0;
 144     Float_t effHighRawEt = 0;
 145     Float_t effLowRawEt = 0;
 146
 147     Float_t nChargedHadronsMeasured = 0.0;
 148     Float_t nChargedHadronsTotal = 0.0;
 149     Float_t nChargedHadronsEtMeasured = 0.0;
 150     Float_t nChargedHadronsEtTotal = 0.0;
 151
 152
 153     for (Int_t iCluster = 0; iCluster < event->GetNumberOfCaloClusters(); iCluster++)
 154     {
 155         AliESDCaloCluster* cluster = event->GetCaloCluster(iCluster);
 156         if (!cluster)
 157         {
 158             AliError(Form("ERROR: Could not get cluster %d", iCluster));
 159             continue;
 160         }
 161         int x = 0;
 162         fCutFlow->Fill(x++);
 163         if(!fSelector->IsDetectorCluster(*cluster)) continue;
 164         fCutFlow->Fill(x++);
 165         if(!fSelector->PassMinEnergyCut(*cluster)) continue;
 166         fCutFlow->Fill(x++);
 167         if (!fSelector->PassDistanceToBadChannelCut(*cluster)) continue;
 168         fCutFlow->Fill(x++);
 169
 170         Float_t pos[3];
 171
 172         cluster->GetPosition(pos);
 173         TVector3 cp(pos);
 174
 175         Bool_t matched = kTRUE;//default to no track matched
 176         Int_t trackMatchedIndex = cluster->GetTrackMatchedIndex();//find the index of the matched track
 177         matched = fSelector->PassTrackMatchingCut(*cluster);
 178         if(!matched){
 179           if(trackMatchedIndex < 0) matched=kTRUE;
 180           AliESDtrack *track = event->GetTrack(trackMatchedIndex);
 181           //if this is a good track, accept track will return true.  The track matched is good, so not track matched is false
 182           matched = !(fEsdtrackCutsTPC->AcceptTrack(track));
 183         }
 184
 185
 186         if (matched)
 187         {
 188
 189             if (cluster->GetNTracksMatched() > 0 && trackMatchedIndex>=0)
 190             {
 191                 AliVTrack *track = event->GetTrack(trackMatchedIndex);
 192                 if (!track) {
 193                     AliError("Error: track does not exist");
 194                 }
 195                 else {
 196                   nChargedHadronsMeasured++;
 197                   nChargedHadronsTotal += 1/fTmCorrections->TrackMatchingEfficiency(track->Pt(),fClusterMult);
 198                   Double_t effCorrEt = CorrectForReconstructionEfficiency(*cluster,fClusterMult);
 199                   nChargedHadronsEtMeasured+= TMath::Sin(cp.Theta())*effCorrEt;
 200                   nChargedHadronsEtTotal+= 1/fTmCorrections->TrackMatchingEfficiency(track->Pt(),fClusterMult) *effCorrEt;
 201                   fHistMatchedTracksEvspTvsMult->Fill(track->P(),TMath::Sin(cp.Theta())*cluster->E(),fClusterMult);
 202                   fHistMatchedTracksEvspTvsMultEffCorr->Fill(track->P(),CorrectForReconstructionEfficiency(*cluster,fClusterMult),fClusterMult);
 203                     const Double_t *pidWeights = track->PID();
 204
 205                     Double_t maxpidweight = 0;
 206                     Int_t maxpid = 0;
 207
 208                     if (pidWeights)
 209                     {
 210                         for (Int_t p =0; p < AliPID::kSPECIES; p++)
 211                         {
 212                             if (pidWeights[p] > maxpidweight)
 213                             {
 214                                 maxpidweight = pidWeights[p];
 215                                 maxpid = p;
 216                             }
 217                         }
 218                         if (fCuts->GetHistMakeTreeDeposit() && fDepositTree)
 219                         {
 220                             fEnergyDeposited = cluster->E();
 221                             fMomentumTPC = track->P();
 222                             fCharge = track->Charge();
 223                             fParticlePid = maxpid;
 224                             fPidProb = maxpidweight;
 225                             AliESDtrack *esdTrack = dynamic_cast<AliESDtrack*>(track);
 226                             if (!esdTrack) {
 227                                 AliError("Error: track does not exist");
 228                             }
 229                             else {
 230                                 if (esdTrack) fTrackPassedCut = fEsdtrackCutsTPC->AcceptTrack(esdTrack);
 231                                 fDepositTree->Fill();
 232                             }
 233                         }
 234
 235                         if (maxpidweight > fPidCut)
 236                         {
 237                             //Float_t dist = TMath::Sqrt(pos[0]*pos[0] + pos[1]*pos[1]);
 238
 239                             //Float_t theta = TMath::ATan(pos[2]/dist)+TMath::Pi()/2;
 240
 241                             //Float_t et = cluster->E() * TMath::Sin(theta);
 242                             if (maxpid == AliPID::kProton)
 243                             {
 244
 245                                 if (track->Charge() == 1)
 246                                 {
 247                                     fHistProtonEnergyDeposit->Fill(cluster->E(), track->E());
 248                                 }
 249                                 else if (track->Charge() == -1)
 250                                 {
 251                                     fHistAntiProtonEnergyDeposit->Fill(cluster->E(), track->E());
 252                                 }
 253                             }
 254                             else if (maxpid == AliPID::kPion)
 255                             {
 256                                 fHistChargedPionEnergyDeposit->Fill(cluster->E(), track->E());
 257                             }
 258                             else if (maxpid == AliPID::kKaon)
 259                             {
 260                                 fHistChargedKaonEnergyDeposit->Fill(cluster->E(), track->E());
 261                             }
 262                             else if (maxpid == AliPID::kMuon)
 263                             {
 264                                 fHistMuonEnergyDeposit->Fill(cluster->E(), track->E());
 265                             }
 266                         }
 267                     }
 268                 }
 269             }
 270             //continue;
 271         } // distance
 272         else{//these are clusters which were not track matched
 273           fCutFlow->Fill(x++);
 274           //std::cout << x++ << std::endl;
 275
 276           //if (cluster->E() >  fSingleCellEnergyCut && cluster->GetNCells() == fCuts->GetCommonSingleCell()) continue;
 277           //if (cluster->E() < fClusterEnergyCut) continue;
 278           cluster->GetPosition(pos);
 279
 280             TVector3 p2(pos);
 281
 282             fClusterPosition->Fill(p2.Phi(), p2.PseudoRapidity());
 283             fClusterEnergy->Fill(cluster->E());
 284             fClusterEt->Fill(TMath::Sin(p2.Theta())*cluster->E());
 285
 286             Double_t effCorrEt = CorrectForReconstructionEfficiency(*cluster,fClusterMult);
 287             fTotNeutralEt += effCorrEt;
 288             nominalRawEt += effCorrEt;
 289             nonlinHighRawEt += effCorrEt*GetCorrectionModification(*cluster,1,0);
 290             nonlinLowRawEt += effCorrEt*GetCorrectionModification(*cluster,-1,0);
 291             effHighRawEt += effCorrEt*GetCorrectionModification(*cluster,0,1);
 292             effLowRawEt += effCorrEt*GetCorrectionModification(*cluster,0,-1);
 293             fNeutralMultiplicity++;
 294         }
 295         fMultiplicity++;
 296     }
 297
 298     fChargedEnergyRemoved = GetChargedContribution(fNeutralMultiplicity);
 299     fNeutralEnergyRemoved = GetNeutralContribution(fNeutralMultiplicity);
 300     fHistChargedEnergyRemoved->Fill(fChargedEnergyRemoved, fNeutralMultiplicity);
 301     fHistNeutralEnergyRemoved->Fill(fNeutralEnergyRemoved, fNeutralMultiplicity);
 302
 303     fGammaEnergyAdded = GetGammaContribution(fNeutralMultiplicity);
 304     fHistGammaEnergyAdded->Fill(fGammaEnergyAdded, fNeutralMultiplicity);
 305
 306     Double_t removedEnergy = GetChargedContribution(fNeutralMultiplicity) + GetNeutralContribution(fNeutralMultiplicity) + GetGammaContribution(fNeutralMultiplicity) + GetSecondaryContribution(fNeutralMultiplicity);
 307     fHistRemovedEnergy->Fill(removedEnergy);
 308
 309     fTotNeutralEt = fGeomCorrection * fEMinCorrection * (fTotNeutralEt - removedEnergy);
 310     fTotEt = fTotChargedEt + fTotNeutralEt;
 311 // Fill the histograms...0
 312     FillHistograms();
 313     //std::cout << "fTotNeutralEt: " << fTotNeutralEt << ", Contribution from non-removed charged: " << GetChargedContribution(fNeutralMultiplicity) << ", neutral: " << GetNeutralContribution(fNeutralMultiplicity) << ", gammas: " << GetGammaContribution(fNeutralMultiplicity) << ", multiplicity: " << fNeutralMultiplicity<< std::endl;
 314     //cout<<"cent "<<cent<<" cluster mult "<<fClusterMult<<" fTotNeutralEt "<<fTotNeutralEt<<" nominalRawEt "<<nominalRawEt<<endl;
 315     fHistNominalRawEt->Fill(nominalRawEt,cent);
 316     fHistNominalNonLinHighEt->Fill(nonlinHighRawEt,cent);
 317     fHistNominalNonLinLowEt->Fill(nonlinLowRawEt,cent);
 318     fHistNominalEffHighEt->Fill(effHighRawEt,cent);
 319     fHistNominalEffLowEt->Fill(effLowRawEt,cent);
 320     fHistFoundHadronsvsCent->Fill(nChargedHadronsMeasured,cent);
 321     fHistNotFoundHadronsvsCent->Fill(nChargedHadronsTotal-nChargedHadronsMeasured,cent);
 322     fHistFoundHadronsEtvsCent->Fill(nChargedHadronsEtMeasured,cent);
 323     fHistNotFoundHadronsEtvsCent->Fill(nChargedHadronsEtTotal-nChargedHadronsEtMeasured,cent);
 324 //     cout<<"Number of hadrons measured:  "<<nChargedHadronsMeasured<<" Estimated total number of hadrons "<<nChargedHadronsTotal<<" ET in track matched hadrons "<<
 325 //       nChargedHadronsEtMeasured;
 326 //     if(nChargedHadronsMeasured>0)cout<<" ("<<nChargedHadronsEtMeasured/nChargedHadronsMeasured<<") ";
 327 //     cout<<" ET in all hadrons ";
 328 //     cout<<nChargedHadronsEtTotal;
 329 //     if(nChargedHadronsTotal>0) cout<<" ("<<nChargedHadronsEtTotal/nChargedHadronsTotal<<") ";
 330 //     cout<<endl;
 331     return 0;
 332 }
 333
 334 bool AliAnalysisEtReconstructed::CheckGoodVertex(AliVParticle* track)
 335 { // check vertex
 336
 337     Float_t bxy = 999.;
 338     Float_t bz = 999.;
 339     if (!track) {
 340         AliError("ERROR: no track");
 341         return kFALSE;
 342     }
 343     AliESDtrack *esdTrack = dynamic_cast<AliESDtrack*>(track);
 344     if (!esdTrack) {
 345         AliError("ERROR: no track");
 346         return kFALSE;
 347     }
 348     esdTrack->GetImpactParametersTPC(bxy,bz);
 349
 350
 351     bool status = (TMath::Abs(track->Xv()) < fCuts->GetReconstructedVertexXCut()) &&
 352                   (TMath::Abs(track->Yv()) < fCuts->GetReconstructedVertexYCut()) &&
 353                   (TMath::Abs(track->Zv()) < fCuts->GetReconstructedVertexZCut()) &&
 354                   (TMath::Abs(bxy) < fCuts->GetReconstructedIPxyCut()) &&
 355                   (TMath::Abs(bz) < fCuts->GetReconstructedIPzCut());
 356
 357     return status;
 358 }
 359
 360 void AliAnalysisEtReconstructed::Init()
 361 { // Init
 362     AliAnalysisEt::Init();
 363     fPidCut = fCuts->GetReconstructedPidCut();
 364     TGeoGlobalMagField::Instance()->SetField(new AliMagF("Maps","Maps", 1., 1., AliMagF::k5kG));
 365     if (!fCorrections) {
 366         cout<<"Warning!  You have not set corrections.  Your code will crash.  You have to set the corrections."<<endl;
 367     }
 368 }
 369
 370 bool AliAnalysisEtReconstructed::TrackHitsCalorimeter(AliVParticle* track, Double_t magField)
 371 { // propagate track to detector radius
 372
 373     if (!track) {
 374         cout<<"Warning: track empty"<<endl;
 375         return kFALSE;
 376     }
 377     AliESDtrack *esdTrack= dynamic_cast<AliESDtrack*>(track);
 378     if (!esdTrack) {
 379         AliError("ERROR: no ESD track");
 380         return kFALSE;
 381     }
 382     // Printf("Propagating track: eta: %f, phi: %f, pt: %f", esdTrack->Eta(), esdTrack->Phi(), esdTrack->Pt());
 383
 384     Bool_t prop = esdTrack->PropagateTo(fDetectorRadius, magField);
 385
 386     // if (prop) Printf("Track propagated, eta: %f, phi: %f, pt: %f", esdTrack->Eta(), esdTrack->Phi(), esdTrack->Pt());
 387     return prop && fSelector->CutGeometricalAcceptance(*esdTrack);
 388 }
 389
 390 void AliAnalysisEtReconstructed::FillOutputList(TList* list)
 391 { // add some extra histograms to the ones from base class
 392     AliAnalysisEt::FillOutputList(list);
 393
 394     list->Add(fHistChargedPionEnergyDeposit);
 395     list->Add(fHistProtonEnergyDeposit);
 396     list->Add(fHistAntiProtonEnergyDeposit);
 397     list->Add(fHistChargedKaonEnergyDeposit);
 398     list->Add(fHistMuonEnergyDeposit);
 399
 400     list->Add(fHistRemovedEnergy);
 401     list->Add(fClusterPosition);
 402     list->Add(fClusterEnergy);
 403     list->Add(fClusterEt);
 404
 405     list->Add(fHistChargedEnergyRemoved);
 406     list->Add(fHistNeutralEnergyRemoved);
 407     list->Add(fHistGammaEnergyAdded);
 408     list->Add(fHistMatchedTracksEvspTvsMult);
 409     list->Add(fHistMatchedTracksEvspTvsMultEffCorr);
 410     list->Add(fHistFoundHadronsvsCent);
 411     list->Add(fHistNotFoundHadronsvsCent);
 412     list->Add(fHistFoundHadronsEtvsCent);
 413     list->Add(fHistNotFoundHadronsEtvsCent);
 414     list->Add(fHistNominalRawEt);
 415     list->Add(fHistNominalNonLinHighEt);
 416     list->Add(fHistNominalNonLinLowEt);
 417     list->Add(fHistNominalEffHighEt);
 418     list->Add(fHistNominalEffLowEt);
 419 }
 420
 421 void AliAnalysisEtReconstructed::CreateHistograms()
 422 { // add some extra histograms to the ones from base class
 423     AliAnalysisEt::CreateHistograms();
 424
 425     Int_t nbinsEt = 1000;
 426     Double_t minEt = 0;
 427     Double_t maxEt = 10;
 428
 429     // possibly change histogram limits
 430 //     if (fCuts) {
 431 //         nbinsEt = fCuts->GetHistNbinsParticleEt();
 432 //         minEt = fCuts->GetHistMinParticleEt();
 433 //         maxEt = fCuts->GetHistMaxParticleEt();
 434 //     }
 435
 436     TString histname;
 437     histname = "fHistChargedPionEnergyDeposit" + fHistogramNameSuffix;
 438     fHistChargedPionEnergyDeposit = new TH2F(histname.Data(), "Energy deposited by #pi^{+/-}", nbinsEt, minEt, maxEt, nbinsEt, minEt, maxEt);
 439     fHistChargedPionEnergyDeposit->SetXTitle("Energy deposited in calorimeter");
 440     fHistChargedPionEnergyDeposit->SetYTitle("Energy of track");
 441
 442     histname = "fHistProtonEnergyDeposit" + fHistogramNameSuffix;
 443     fHistProtonEnergyDeposit = new TH2F(histname.Data(), "Energy deposited by protons", nbinsEt, minEt, maxEt, nbinsEt, minEt, maxEt);
 444     fHistProtonEnergyDeposit->SetXTitle("Energy deposited in calorimeter");
 445     fHistProtonEnergyDeposit->SetYTitle("Energy of track");
 446
 447     histname = "fHistAntiProtonEnergyDeposit" + fHistogramNameSuffix;
 448     fHistAntiProtonEnergyDeposit = new TH2F(histname.Data(), "Energy deposited by anti-protons", nbinsEt, minEt, maxEt, nbinsEt, minEt, maxEt);
 449     fHistAntiProtonEnergyDeposit->SetXTitle("Energy deposited in calorimeter");
 450     fHistAntiProtonEnergyDeposit->SetYTitle("Energy of track");
 451
 452     histname = "fHistChargedKaonEnergyDeposit" + fHistogramNameSuffix;
 453     fHistChargedKaonEnergyDeposit = new TH2F(histname.Data(), "Energy deposited by K^{+/-}", nbinsEt, minEt, maxEt, nbinsEt, minEt, maxEt);
 454     fHistChargedKaonEnergyDeposit->SetXTitle("Energy deposited in calorimeter");
 455     fHistChargedKaonEnergyDeposit->SetYTitle("Energy of track");
 456
 457     histname = "fHistMuonEnergyDeposit" + fHistogramNameSuffix;
 458     fHistMuonEnergyDeposit = new TH2F(histname.Data(), "Energy deposited by #mu^{+/-}", nbinsEt, minEt, maxEt, nbinsEt, minEt, maxEt);
 459     fHistMuonEnergyDeposit->SetXTitle("Energy deposited in calorimeter");
 460     fHistMuonEnergyDeposit->SetYTitle("Energy of track");
 461
 462     histname = "fHistRemovedEnergy" + fHistogramNameSuffix;
 463     fHistRemovedEnergy = new TH1F(histname.Data(), histname.Data(), 1000, 0, 20);
 464     //fHistMuonEnergyDeposit->SetXTitle("Energy deposited in calorimeter");
 465     //fHistMuonEnergyDeposit->SetYTitle("Energy of track");
 466
 467     histname = "fClusterPosition" + fHistogramNameSuffix;
 468     fClusterPosition = new TH2D(histname.Data(), "Position of accepted neutral clusters",1000, -2.0, -.5, 1000, -.13 , 0.13);
 469     fClusterPosition->SetXTitle("Energy deposited in calorimeter");
 470     fClusterPosition->SetYTitle("Energy of track");
 471
 472     histname = "fClusterEnergy" + fHistogramNameSuffix;
 473     fClusterEnergy = new TH1F(histname.Data(), histname.Data(), 100, 0, 5);
 474     fClusterEnergy->SetXTitle("Number of clusters");
 475     fClusterEnergy->SetYTitle("Energy of cluster");
 476
 477     histname = "fClusterEt" + fHistogramNameSuffix;
 478     fClusterEt = new TH1F(histname.Data(), histname.Data(), 100, 0, 5);
 479     fClusterEt->SetXTitle("Number of clusters");
 480     fClusterEt->SetYTitle("E_{T} of cluster");
 481
 482     histname = "fHistChargedEnergyRemoved" + fHistogramNameSuffix;
 483     fHistChargedEnergyRemoved = new TH2D(histname.Data(), histname.Data(), 1000, .0, 30, 100, -0.5 , 99.5);
 484
 485     histname = "fHistNeutralEnergyRemoved" + fHistogramNameSuffix;
 486     fHistNeutralEnergyRemoved = new TH2D(histname.Data(), histname.Data(), 1000, .0, 30, 100, -0.5 , 99.5);
 487
 488     histname = "fHistGammaEnergyAdded" + fHistogramNameSuffix;
 489     fHistGammaEnergyAdded = new TH2D(histname.Data(), histname.Data(), 1000, .0, 30, 100, -0.5 , 99.5);
 490
 491     fHistMatchedTracksEvspTvsMult = new TH3F("fHistMatchedTracksEvspTvsMult", "fHistMatchedTracksEvspTvsMult",100, 0, 3,100,0,3,10,0,100);
 492     fHistMatchedTracksEvspTvsMultEffCorr = new TH3F("fHistMatchedTracksEvspTvsMultEffCorr", "fHistMatchedTracksEvspTvsMultEffCorr",100, 0, 3,100,0,3,10,0,100);
 493     fHistFoundHadronsvsCent = new TH2F("fHistFoundHadronsvsCent","fHistFoundHadronsvsCent",100,0,100,20,0,20);
 494     fHistNotFoundHadronsvsCent = new TH2F("fHistNotFoundHadronsvsCent","fHistNotFoundHadronsvsCent",100,0,100,20,0,20);
 495     fHistFoundHadronsEtvsCent = new TH2F("fHistFoundHadronsEtvsCent","fHistFoundHadronsEtvsCent",100,0,200,20,0,20);
 496     fHistNotFoundHadronsEtvsCent = new TH2F("fHistNotFoundHadronsEtvsCent","fHistNotFoundHadronsEtvsCent",100,0,200,20,0,20);
 497
 498     maxEt = 100;
 499     histname = "fHistNominalRawEt" + fHistogramNameSuffix;
 500     fHistNominalRawEt = new TH2D(histname.Data(), histname.Data(),nbinsEt,minEt,maxEt,20,-0.5,19.5);
 501     histname = "fHistNominalNonLinHighEt" + fHistogramNameSuffix;
 502     fHistNominalNonLinHighEt = new TH2D(histname.Data(), histname.Data(),nbinsEt,minEt,maxEt,20,-0.5,19.5);
 503     histname = "fHistNominalNonLinLowEt" + fHistogramNameSuffix;
 504     fHistNominalNonLinLowEt = new TH2D(histname.Data(), histname.Data(),nbinsEt,minEt,maxEt,20,-0.5,19.5);
 505     histname = "fHistNominalEffHighEt" + fHistogramNameSuffix;
 506     fHistNominalEffHighEt = new TH2D(histname.Data(), histname.Data(),nbinsEt,minEt,maxEt,20,-0.5,19.5);
 507     histname = "fHistNominalEffLowEt" + fHistogramNameSuffix;
 508     fHistNominalEffLowEt = new TH2D(histname.Data(), histname.Data(),nbinsEt,minEt,maxEt,20,-0.5,19.5);
 509
 510 }
 511 Double_t AliAnalysisEtReconstructed::ApplyModifiedCorrections(const AliESDCaloCluster& cluster,Int_t nonLinCorr, Int_t effCorr, Int_t mult)
 512 {
 513   Float_t pos[3];
 514   cluster.GetPosition(pos);
 515   TVector3 cp(pos);
 516   Double_t corrEnergy = fReCorrections->CorrectedEnergy(cluster.E(),mult);
 517
 518   Double_t factorNonLin = GetCorrectionModification(cluster, nonLinCorr,effCorr);
 519
 520   //std::cout << "Original energy: " << cluster.E() << ", corrected energy: " << corrEnergy << std::endl;
 521   return TMath::Sin(cp.Theta())*corrEnergy*factorNonLin;
 522 }
 523
 524 Double_t AliAnalysisEtReconstructed::GetCorrectionModification(const AliESDCaloCluster& cluster,Int_t nonLinCorr, Int_t effCorr){//nonLinCorr 0 = nominal 1 = high -1 = low, effCorr  0 = nominal 1 = high -1 = low
 525   if(nonLinCorr==0){
 526     cout<<"Warning:  This function should not get called!"<<endl;//this statement is basically here to avoid a compilation warning
 527   }
 528   if(effCorr==0){
 529     cout<<"Warning:  This function should not get called!"<<endl;//this statement is basically here to avoid a compilation warning
 530   }
 531   return cluster.E();
 532 }