PWGGA/PHOSTasks/PHOS_EpRatio/AliAnalysisTaskEpRatio.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 #include "TObjArray.h"
  17 #include "TF1.h"
  18 #include "TH1F.h"
  19 #include "TH2F.h"
  20 #include "TH2I.h"
  21 #include "TH3F.h"
  22 #include "TParticle.h"
  23 #include "TList.h"
  24 #include "THashList.h"
  25
  26 #include "AliAnalysisTaskSE.h"
  27 #include "AliAnalysisTaskEpRatio.h"
  28 #include "AliAnalysisManager.h"
  29 #include "AliInputEventHandler.h"
  30 #include "AliCaloPhoton.h"
  31 #include "AliPHOSGeometry.h"
  32 #include "AliVEvent.h"
  33 #include "AliAODEvent.h"
  34 #include "AliVCaloCells.h"
  35 #include "AliVCluster.h"
  36 #include "AliVTrack.h"
  37 #include "AliLog.h"
  38 #include "AliPID.h"
  39 #include "AliPIDResponse.h"
  40
  41 // *********************** New
  42 #include "AliVertex.h"
  43 #include "AliVVertex.h"
  44 #include "AliVertexerTracks.h"
  45 #include "AliExternalTrackParam.h"
  46 #include "AliTrackReference.h"
  47
  48 #include <AliCDBManager.h>
  49 #include "AliMagF.h"
  50 #include "TGeoGlobalMagField.h"
  51 #include "AliPHOSCalibData.h"
  52 #include "AliOADBContainer.h"
  53
  54 // Analysis task to fill histograms with E/p ratios for electrons, positrons, muons
  55 // and hadrons (pions, kaons..), where E is an energy of PHOS cluster
  56 // and p is the momentum of it's matched track.
  57
  58 // Authors: Boris Polishchuk,Tsubasa Okubo
  59 // Date   : 04.07.2013
  60
  61 ClassImp(AliAnalysisTaskEpRatio)
  62 //________________________________________________________________________
  63 AliAnalysisTaskEpRatio::AliAnalysisTaskEpRatio(const char *name) :  AliAnalysisTaskSE(name),
  64   fRunNumber(-999),
  65   fOutputContainer(0),
  66   fPHOSGeo(0),
  67   fPIDResponse(0) // Yuri
  68 // Output slots #0 write into a TH1 container
  69 {
  70   DefineOutput(1,TList::Class());
  71
  72 }
  73 //________________________________________________________________________
  74 void AliAnalysisTaskEpRatio::UserCreateOutputObjects()
  75 {
  76   // Create histograms
  77   // Called once
  78
  79   if(fOutputContainer != NULL){
  80     delete fOutputContainer;
  81   }
  82   fOutputContainer = new THashList();
  83   fOutputContainer->SetOwner();
  84
  85   fOutputContainer->Add( new TH1F("h_CellMultEvent","PHOS Cell Multiplicity per Event",200,0,200) );
  86
  87   fOutputContainer->Add( new TH2F("hEp","E/p All charged particles (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
  88   fOutputContainer->Add( new TH2F("hEp_ele","E/p electron & positron (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
  89   fOutputContainer->Add( new TH2F("hEp_ele2","E/p electron & positron (E/p VS. 1/p)",1000,0,2.,500,0,5) );
  90   fOutputContainer->Add( new TH2F("hEp_other","E/p other (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
  91
  92   fOutputContainer->Add( new TH2F("hEp_electron","E/p Electron only (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
  93   fOutputContainer->Add( new TH2F("hEp_positron","E/p Positron only (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
  94
  95   fOutputContainer->Add( new TH2F("hEp_mod1","E/p All charged particles in Module1 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
  96   fOutputContainer->Add( new TH2F("hEp_mod2","E/p All charged particles in Module2 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
  97   fOutputContainer->Add( new TH2F("hEp_mod3","E/p All charged particles in Module3 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
  98
  99   fOutputContainer->Add( new TH2F("hEp_ele_mod1","E/p electron & positron in Module1 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
 100   fOutputContainer->Add( new TH2F("hEp_ele_mod2","E/p electron & positron in Module2 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
 101   fOutputContainer->Add( new TH2F("hEp_ele_mod3","E/p electron & positron in Module3 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
 102
 103   fOutputContainer->Add( new TH2F("hEp_other_mod1","E/p other in Module1 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
 104   fOutputContainer->Add( new TH2F("hEp_other_mod2","E/p other in Module2 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
 105   fOutputContainer->Add( new TH2F("hEp_other_mod3","E/p other in Module3 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
 106
 107   fOutputContainer->Add( new TH2F("h_dedx_PHOS","dE/dx in PHOS via TPC",1000,0,10,2200,0,220.) );
 108   fOutputContainer->Add( new TH2F("h_dedx_ele_PHOS","dE/dx electron in PHOS via TPC",1000,0,10,2200,0,220.) );
 109   fOutputContainer->Add( new TH2F("h_dedx_other_PHOS","dE/dx other in PHOS via TPC",1000,0,10,2200,0,220.) );
 110   fOutputContainer->Add( new TH2F("h_dedx_pion_PHOS"  ,"dE/dx pion in PHOS via TPC",1000,0,10,2200,0,220.) );
 111   fOutputContainer->Add( new TH2F("h_dedx_proton_PHOS","dE/dx proton in PHOS via TPC",1000,0,10,2200,0,220.) );
 112   fOutputContainer->Add( new TH2F("h_dedx_kaon_PHOS"  ,"dE/dx kaon in PHOS via TPC",1000,0,10,2200,0,220.) );
 113   fOutputContainer->Add( new TH2F("h_dedx_muon_PHOS"  ,"dE/dx muon in PHOS via TPC",1000,0,10,2200,0,220.) );
 114
 115   fOutputContainer->Add( new TH2F("h_ClusNXZM1","Cluster (X,Z) Module1",   64,0.5,64.5, 56,0.5,56.5) );
 116   fOutputContainer->Add( new TH2F("h_ClusNXZM2","Cluster (X,Z) Module2",   64,0.5,64.5, 56,0.5,56.5) );
 117   fOutputContainer->Add( new TH2F("h_ClusNXZM3","Cluster (X,Z) Module3",   64,0.5,64.5, 56,0.5,56.5) );
 118   fOutputContainer->Add( new TH2F("h_ClusNXZM1_ele","Electron Cluster (X,Z) Module1", 64,0.5,64.5, 56,0.5,56.5) );
 119   fOutputContainer->Add( new TH2F("h_ClusNXZM2_ele","Electron Cluster (X,Z) Module2", 64,0.5,64.5, 56,0.5,56.5) );
 120   fOutputContainer->Add( new TH2F("h_ClusNXZM3_ele","Electron Cluster (X,Z) Module3", 64,0.5,64.5, 56,0.5,56.5) );
 121
 122   fOutputContainer->Add( new TH2F("h_EClus_NCell","Energy VS. NCell in PHOS",1000,0,10,30,0,30.) );
 123   fOutputContainer->Add( new TH1F("h_ECluster","Energy of Cluster in PHOS",1000,0,10.) );
 124   fOutputContainer->Add( new TH1F("h_NCells","Number of Cells per Cluster in PHOS",30,0,30.) );
 125
 126   // =====================================================================================================
 127   fOutputContainer->Add( new TH2F("hEp_had","E/p Hadrons (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
 128   fOutputContainer->Add( new TH2F("hEp_had2","E/p Hadrons (E/p VS. 1/p)",1000,0,2.,500,0,5) );
 129   fOutputContainer->Add( new TH2F("h_dedx_had_PHOS","dE/dx Hadron in PHOS via TPC",1000,0,10,2200,0,220.) );
 130   fOutputContainer->Add( new TH2F("h_ClusNXZM1_had","Hadron Cluster (X,Z) Module1", 64,0.5,64.5, 56,0.5,56.5) );
 131   fOutputContainer->Add( new TH2F("h_ClusNXZM2_had","Hadron Cluster (X,Z) Module2", 64,0.5,64.5, 56,0.5,56.5) );
 132   fOutputContainer->Add( new TH2F("h_ClusNXZM3_had","Hadron Cluster (X,Z) Module3", 64,0.5,64.5, 56,0.5,56.5) );
 133   fOutputContainer->Add( new TH2F("hEp_had_mod1","E/p Hadron in Module1 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
 134   fOutputContainer->Add( new TH2F("hEp_had_mod2","E/p Hadron in Module2 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
 135   fOutputContainer->Add( new TH2F("hEp_had_mod3","E/p Hadron in Module3 (E/p VS. p_{T})",1000,0,2,1000,0,10.) );
 136
 137   // =====================================================================================================
 138   fOutputContainer->Add( new TH1F("h_NClusEvent","Number PHOS Clusters per Event in All Modules",31,-0.5,30.5) );
 139   fOutputContainer->Add( new TH1F("h_NClusEventmod1","Number PHOS Clusters per Event in Module1",21,-0.5,20.5) );
 140   fOutputContainer->Add( new TH1F("h_NClusEventmod2","Number PHOS Clusters per Event in Module2",21,-0.5,20.5) );
 141   fOutputContainer->Add( new TH1F("h_NClusEventmod3","Number PHOS Clusters per Event in Module3",21,-0.5,20.5) );
 142
 143   fOutputContainer->Add( new TH1F("h_EClusEvent","PHOS Cluster Energy per Event in All Modules",2010,-0.5,20.5) );
 144   fOutputContainer->Add( new TH1F("h_EClusEventmod1","PHOS Cluster Energy per Event in Module1",2010,-0.5,20.5) );
 145   fOutputContainer->Add( new TH1F("h_EClusEventmod2","PHOS Cluster Energy per Event in Module2",2010,-0.5,20.5) );
 146   fOutputContainer->Add( new TH1F("h_EClusEventmod3","PHOS Cluster Energy per Event in Module3",2010,-0.5,20.5) );
 147
 148   fOutputContainer->Add( new TH1F("h_NTrackEvent","Number of TPC Tracks per Event",1001,-0.5,1000.5) );
 149   fOutputContainer->Add( new TH1F("h_PTrackEvent","TPC Momentum per Event",2010,-0.5,20.5) );
 150   fOutputContainer->Add( new TH1F("h_NClusTPCEvent","Number of TPC Clusters per Event",301,-0.5,300.5) );
 151
 152   // ======================================================================================================
 153
 154
 155   AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 156   AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
 157
 158   fPIDResponse = inputHandler->GetPIDResponse();
 159   if(!fPIDResponse) AliFatal(" !!! FATAL!! No PIDResponse!!");
 160
 161   PostData(1, fOutputContainer);
 162
 163 }
 164
 165 //________________________________________________________________________
 166 void AliAnalysisTaskEpRatio::UserExec(Option_t *)
 167 {
 168   //Main function that does all the job.
 169
 170   AliVEvent* event = InputEvent();
 171   fRunNumber = event->GetRunNumber();
 172
 173   SetGeometry();
 174
 175   Double_t v[3]={0,0,0}; //vertex;
 176   const AliVVertex *primvertex = event->GetPrimaryVertex();
 177
 178   primvertex->GetXYZ(v);
 179   TVector3 vtx(v);
 180
 181   if(fabs(v[2])>10.) return;
 182
 183   const Int_t kNumberOfTracks = event->GetNumberOfTracks();
 184   const Int_t kNumberOfCaloClusters = event->GetNumberOfCaloClusters();
 185
 186   AliVCaloCells *cells = event->GetPHOSCells();
 187   Int_t multCells = cells->GetNumberOfCells();
 188   FillHistogram("h_CellMultEvent",multCells);
 189
 190   TLorentzVector pc1; //4-momentum of PHOS cluster 1
 191   Int_t inPHOS=0;
 192
 193   // ================
 194   // === TPC dedx ===
 195   // ================
 196   FillHistogram("h_NTrackEvent",kNumberOfTracks);
 197
 198   for(Int_t iTrack=0; iTrack<kNumberOfTracks; iTrack++){
 199     AliVTrack* esdTrack = dynamic_cast<AliVTrack*> (event->GetTrack(iTrack));
 200     FillHistogram("h_PTrackEvent",esdTrack->P() );
 201     FillHistogram("h_NClusTPCEvent",esdTrack->GetTPCNcls() );
 202   }
 203
 204   // ================
 205   // === E/p ========
 206   // ================
 207   Float_t position[3];
 208   Int_t mod1, relId[4], cellAbsId, cellX, cellZ;
 209   Int_t nPHOSCluster = 0, inMod1 = 0, inMod2 = 0, inMod3 = 0;
 210
 211   for(Int_t iClu=0; iClu<kNumberOfCaloClusters; iClu++){
 212
 213     AliVCluster *c1 = event->GetCaloCluster(iClu);
 214
 215     if( !c1->IsPHOS() ) continue;
 216     //if( c1->E()<0.3 ) continue;
 217     FillHistogram("h_EClus_NCell",c1->E(),c1->GetNCells());
 218     FillHistogram("h_ECluster",c1->E());
 219     FillHistogram("h_NCells",c1->GetNCells());
 220     if( c1->GetNCells()<3 ) continue;
 221
 222     c1->GetPosition(position);
 223     TVector3 global1(position);
 224     fPHOSGeo->GlobalPos2RelId(global1,relId);
 225     mod1  = relId[0];
 226     cellX = relId[2];
 227     cellZ = relId[3];
 228
 229     cellAbsId = c1->GetCellAbsId(0);
 230     fPHOSGeo->AbsToRelNumbering(cellAbsId,relId);
 231     mod1 = relId[0];
 232
 233     FillHistogram("h_EClusEvent",c1->E());
 234
 235     if( mod1 == 1 ){
 236       inMod1++; nPHOSCluster++;
 237       FillHistogram("h_EClusEventmod1",c1->E());
 238     }
 239     else if( mod1 == 2 ){
 240       inMod2++; nPHOSCluster++;
 241       FillHistogram("h_EClusEventmod2",c1->E());
 242     }
 243     else if( mod1 == 3 ){
 244       inMod3++; nPHOSCluster++;
 245       FillHistogram("h_EClusEventmod3",c1->E());
 246     }
 247
 248     Int_t nMatched = c1->GetNTracksMatched();
 249
 250     Double_t Dx = c1->GetTrackDx();
 251     Double_t Dz = c1->GetTrackDz();
 252     Double_t r = sqrt(Dx*Dx+Dz*Dz); // unit is [cm]
 253
 254     if(r > 5.) continue;
 255     if(nMatched<=0) continue;
 256
 257     AliVTrack* esdTrack = dynamic_cast<AliVTrack*> (c1->GetTrackMatched(0));
 258     if( !(TMath::Abs(esdTrack->Eta())< 0.8) ) continue;
 259
 260     Short_t charge   = esdTrack->Charge();
 261     if(fPIDResponse) {
 262
 263       //if(c1->E()<0.3) continue;  // MIP & noise cut
 264
 265       // =============================================================================
 266       FillHistogram("hEp",c1->E()/esdTrack->P(),esdTrack->Pt());
 267       FillHistogram("h_dedx_PHOS",esdTrack->P(),esdTrack->GetTPCsignal());
 268
 269       if( mod1 == 1 ){
 270         FillHistogram("hEp_mod1",c1->E()/esdTrack->P(),esdTrack->Pt());
 271         FillHistogram("h_ClusNXZM1",cellX,cellZ,1.);
 272       }
 273       else if( mod1 == 2 ){
 274         FillHistogram("hEp_mod2",c1->E()/esdTrack->P(),esdTrack->Pt());
 275         FillHistogram("h_ClusNXZM2",cellX,cellZ,1.);
 276       }
 277       else if( mod1 == 3 ){
 278         FillHistogram("hEp_mod3",c1->E()/esdTrack->P(),esdTrack->Pt());
 279         FillHistogram("h_ClusNXZM3",cellX,cellZ,1.);
 280       }
 281
 282
 283       Bool_t pidPion=kFALSE, pidKaon=kFALSE, pidProton=kFALSE, pidElectron=kFALSE, pidMuon=kFALSE, pidHadron=kFALSE;
 284       Double_t nsigmaProton   = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(esdTrack, AliPID::kProton)); // esdTrack << inEvHMain
 285       Double_t nsigmaKaon     = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(esdTrack, AliPID::kKaon));
 286       Double_t nsigmaPion     = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(esdTrack, AliPID::kPion));
 287       //Double_t nsigmaElectron = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(esdTrack, AliPID::kElectron));
 288       Double_t nsigmaElectron = fPIDResponse->NumberOfSigmasTPC(esdTrack, AliPID::kElectron);
 289       Double_t nsigmaMuon     = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(esdTrack, AliPID::kMuon));
 290
 291       if((nsigmaKaon    <nsigmaPion) && (nsigmaKaon    <nsigmaProton) && (nsigmaKaon    <nsigmaElectron)
 292          && (nsigmaKaon     <nsigmaMuon    ) && (nsigmaKaon    <3)) pidKaon     = kTRUE;
 293       if((nsigmaPion    <nsigmaKaon) && (nsigmaPion    <nsigmaProton) && (nsigmaPion    <nsigmaElectron)
 294          && (nsigmaPion     <nsigmaMuon    ) && (nsigmaPion    <3)) pidPion     = kTRUE;
 295       if((nsigmaProton  <nsigmaKaon) && (nsigmaProton  <nsigmaPion  ) && (nsigmaProton  <nsigmaElectron)
 296          && (nsigmaProton   <nsigmaMuon    ) && (nsigmaProton  <3)) pidProton   = kTRUE;
 297       if((nsigmaMuon    <nsigmaKaon) && (nsigmaMuon    <nsigmaPion  ) && (nsigmaMuon    <nsigmaProton  )
 298          && (nsigmaMuon     <nsigmaElectron) && (nsigmaMuon    <3)) pidMuon     = kTRUE;
 299       if((nsigmaElectron<nsigmaKaon) && (nsigmaElectron<nsigmaPion  ) && (nsigmaElectron<nsigmaProton  )
 300          && (nsigmaElectron <nsigmaMuon    ) && (nsigmaElectron<3) && (nsigmaElectron>-2)) pidElectron = kTRUE;
 301       if( (nsigmaElectron<-3) ) pidHadron = kTRUE;
 302
 303       if (pidPion){
 304         FillHistogram("h_dedx_pion_PHOS",esdTrack->P(),esdTrack->GetTPCsignal());
 305       }
 306       if (pidKaon){
 307         FillHistogram("h_dedx_kaon_PHOS",esdTrack->P(),esdTrack->GetTPCsignal());
 308       }
 309       if (pidProton){
 310         FillHistogram("h_dedx_proton_PHOS",esdTrack->P(),esdTrack->GetTPCsignal());
 311       }
 312       if (pidMuon){
 313         FillHistogram("h_dedx_muon_PHOS",esdTrack->P(),esdTrack->GetTPCsignal());
 314       }
 315
 316       if (pidElectron){
 317         FillHistogram("hEp_ele",c1->E()/esdTrack->P(),esdTrack->Pt());
 318         FillHistogram("hEp_ele2",c1->E()/esdTrack->P(),1/esdTrack->P());
 319         FillHistogram("h_dedx_ele_PHOS",esdTrack->P(),esdTrack->GetTPCsignal());
 320
 321         // =============================================================================
 322
 323         if (charge>0) FillHistogram("hEp_positron",c1->E()/esdTrack->P(),esdTrack->Pt());
 324         if (charge<0) FillHistogram("hEp_electron",c1->E()/esdTrack->P(),esdTrack->Pt());
 325
 326         if( mod1 == 1 ){
 327           FillHistogram("h_ClusNXZM1_ele",cellX,cellZ,1.);
 328           FillHistogram("hEp_ele_mod1",c1->E()/esdTrack->P(),esdTrack->Pt());
 329         }
 330         else if( mod1 == 2 ){
 331           FillHistogram("h_ClusNXZM2_ele",cellX,cellZ,1.);
 332           FillHistogram("hEp_ele_mod2",c1->E()/esdTrack->P(),esdTrack->Pt());
 333         }
 334         else if( mod1 == 3 ){
 335           FillHistogram("h_ClusNXZM3_ele",cellX,cellZ,1.);
 336           FillHistogram("hEp_ele_mod3",c1->E()/esdTrack->P(),esdTrack->Pt());
 337         }
 338       }
 339       else{
 340         FillHistogram("hEp_other",c1->E()/esdTrack->P(),esdTrack->Pt());
 341         FillHistogram("h_dedx_other_PHOS",esdTrack->P(),esdTrack->GetTPCsignal());
 342
 343         if( mod1 == 1 ){
 344           FillHistogram("hEp_other_mod1",c1->E()/esdTrack->P(),esdTrack->Pt());
 345         }
 346         else if( mod1 == 2 ){
 347           FillHistogram("hEp_other_mod2",c1->E()/esdTrack->P(),esdTrack->Pt());
 348         }
 349         else if( mod1 == 3 ){
 350           FillHistogram("hEp_other_mod3",c1->E()/esdTrack->P(),esdTrack->Pt());
 351         }
 352
 353       }
 354
 355       if (pidHadron){
 356         FillHistogram("hEp_had",c1->E()/esdTrack->P(),esdTrack->Pt());
 357         FillHistogram("hEp_had2",c1->E()/esdTrack->P(),1/esdTrack->P());
 358         FillHistogram("h_dedx_had_PHOS",esdTrack->P(),esdTrack->GetTPCsignal());
 359
 360         // =============================================================================
 361
 362         if( mod1 == 1 ){
 363           FillHistogram("h_ClusNXZM1_had",cellX,cellZ,1.);
 364           FillHistogram("hEp_had_mod1",c1->E()/esdTrack->P(),esdTrack->Pt());
 365         }
 366         else if( mod1 == 2 ){
 367           FillHistogram("h_ClusNXZM2_had",cellX,cellZ,1.);
 368           FillHistogram("hEp_had_mod2",c1->E()/esdTrack->P(),esdTrack->Pt());
 369         }
 370         else if( mod1 == 3 ){
 371           FillHistogram("h_ClusNXZM3_had",cellX,cellZ,1.);
 372           FillHistogram("hEp_had_mod3",c1->E()/esdTrack->P(),esdTrack->Pt());
 373         }
 374       }
 375
 376     }
 377
 378     // =============================================================================================================================
 379     inPHOS++;
 380   }
 381
 382   FillHistogram("h_NClusEvent",nPHOSCluster);
 383   FillHistogram("h_NClusEventmod1",inMod1);
 384   FillHistogram("h_NClusEventmod2",inMod2);
 385   FillHistogram("h_NClusEventmod3",inMod3);
 386
 387   nPHOSCluster = 0;
 388   inMod1 = 0;
 389   inMod2 = 0;
 390   inMod3 = 0;
 391
 392   //delete caloClustersArr;
 393   PostData(1,fOutputContainer);
 394   return ;
 395 }
 396
 397 //_____________________________________________________________________________
 398 void AliAnalysisTaskEpRatio::FillHistogram(const char * key,Double_t x)const{
 399   //FillHistogram
 400   TH1I * tmpI = dynamic_cast<TH1I*>(fOutputContainer->FindObject(key)) ;
 401   if(tmpI){
 402     tmpI->Fill(x) ;
 403     return ;
 404   }
 405   TH1F * tmpF = dynamic_cast<TH1F*>(fOutputContainer->FindObject(key)) ;
 406   if(tmpF){
 407     tmpF->Fill(x) ;
 408     return ;
 409   }
 410   TH1D * tmpD = dynamic_cast<TH1D*>(fOutputContainer->FindObject(key)) ;
 411   if(tmpD){
 412     tmpD->Fill(x) ;
 413     return ;
 414   }
 415   AliInfo(Form("can not find histogram <%s> ",key)) ;
 416 }
 417 //_____________________________________________________________________________
 418 void AliAnalysisTaskEpRatio::FillHistogram(const char * key,Double_t x,Double_t y)const{
 419   //FillHistogram
 420   TObject * tmp = fOutputContainer->FindObject(key) ;
 421   if(!tmp){
 422     AliInfo(Form("can not find histogram <%s> ",key)) ;
 423     return ;
 424   }
 425   if(tmp->IsA() == TClass::GetClass("TH1F")){
 426     ((TH1F*)tmp)->Fill(x,y) ;
 427     return ;
 428   }
 429   if(tmp->IsA() == TClass::GetClass("TH2F")){
 430     ((TH2F*)tmp)->Fill(x,y) ;
 431     return ;
 432   }
 433   AliError(Form("Calling FillHistogram with 2 parameters for histo <%s> of type %s",key,tmp->IsA()->GetName())) ;
 434 }
 435
 436 //_____________________________________________________________________________
 437 void AliAnalysisTaskEpRatio::SetGeometry()
 438 {
 439   //Init geometry.
 440
 441   if(!fPHOSGeo){
 442
 443     AliOADBContainer geomContainer("phosGeo");
 444     geomContainer.InitFromFile("$ALICE_ROOT/OADB/PHOS/PHOSGeometry.root","PHOSRotationMatrixes");
 445     TObjArray *matrixes = (TObjArray*)geomContainer.GetObject(fRunNumber,"PHOSRotationMatrixes");
 446     fPHOSGeo =  AliPHOSGeometry::GetInstance("IHEP") ;
 447     for(Int_t mod=0; mod<5; mod++) {
 448       if(!matrixes->At(mod)) {
 449         if( fDebug )
 450           AliInfo(Form("No PHOS Matrix for mod:%d, geo=%p\n", mod, fPHOSGeo));
 451         continue;
 452       }
 453       else {
 454         fPHOSGeo->SetMisalMatrix(((TGeoHMatrix*)matrixes->At(mod)),mod) ;
 455         if( fDebug >1 )
 456           AliInfo(Form("Adding PHOS Matrix for mod:%d, geo=%p\n", mod, fPHOSGeo));
 457       }
 458     }
 459   }
 460
 461 }