PWG4/PartCorrDep/AliAnaShowerParameter.cxx

   1 /**************************************************************************\r
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
   3  *                                                                        *\r
   4  * Author: The ALICE Off-line Project.                                    *\r
   5  * Contributors are mentioned in the code where appropriate.              *\r
   6  *                                                                        *\r
   7  * Permission to use, copy, modify and distribute this software and its   *\r
   8  * documentation strictly for non-commercial purposes hereby granted      *\r
   9  * without fee, provided that the above copyright notice appears in all   *\r
  10  * copies and that both the copyright notice and this permission notice   *\r
  11  * appear in the supporting documentation. The authors make no claims     *\r
  12  * about the suitability of this software for any purpose. It is          *\r
  13  * provided "as is" without express or implied warranty.                  *\r
  14  **************************************************************************/\r
  15 /* $Id: AliAnaPhoton.cxx 28688 2008-09-11 15:04:07Z gconesab $ */\r
  16 \r
  17 //_________________________________________________________________________\r
  18 //\r
  19 // Class cloned from AliAnaPhoton, main aim is shower shape studies\r
  20 // \r
  21 // \r
  22 //\r
  23 //-- Author: Jocelyn Mlynarz (WSU) and Gustavo Conesa (LPSC-CNRS)\r
  24 //////////////////////////////////////////////////////////////////////////////\r
  25 \r
  26 \r
  27 // --- ROOT system --- \r
  28 #include <TH3F.h>\r
  29 #include <TH2F.h>\r
  30 #include <TClonesArray.h>\r
  31 //#include <TObjString.h>\r
  32 #include <Riostream.h>\r
  33 #include "TParticle.h"\r
  34 //#include <fstream>\r
  35 \r
  36 // --- Analysis system --- \r
  37 #include "AliAnaShowerParameter.h" \r
  38 #include "AliCaloTrackReader.h"\r
  39 #include "AliStack.h"\r
  40 #include "AliCaloPID.h"\r
  41 #include "AliMCAnalysisUtils.h"\r
  42 #include "AliFiducialCut.h"\r
  43 #include "AliAODCaloCluster.h"\r
  44 #include "AliAODMCParticle.h"\r
  45 #include "AliAnalysisManager.h"\r
  46 #include "AliAnalysisTaskParticleCorrelation.h"\r
  47 #include "AliEMCALGeoUtils.h"\r
  48 #include "AliAODEvent.h"\r
  49 \r
  50 \r
  51 ClassImp(AliAnaShowerParameter)\r
  52 \r
  53 //____________________________________________________________________________\r
  54 AliAnaShowerParameter::AliAnaShowerParameter() : \r
  55 AliAnaPartCorrBaseClass(), fCalorimeter(""), fNCellsCutMin(0),\r
  56 fNCellsCutMax(0), fLambdaCut(0), fTimeCutMin(-1), fTimeCutMax(9999999),\r
  57 fhNClusters(0), fhNCellCluster(0), fhEtaPhiPtCluster(0), \r
  58 fhLambdaPtCluster(0), \r
  59 \r
  60 //MC\r
  61 \r
  62 fhLambdaPtPhoton(0), fhLambdaPtPi0(0), fhLambdaPtPion(0), fhPtTruthPi0(0)\r
  63 \r
  64 {\r
  65   //default ctor\r
  66   \r
  67   //Initialize parameters\r
  68   InitParameters();\r
  69   \r
  70 }\r
  71 \r
  72 //____________________________________________________________________________\r
  73 AliAnaShowerParameter::~AliAnaShowerParameter() \r
  74 {\r
  75   //dtor\r
  76   \r
  77 }\r
  78 \r
  79 //________________________________________________________________________\r
  80 TObjString *  AliAnaShowerParameter::GetAnalysisCuts()\r
  81 {       \r
  82   //Save parameters used for analysis\r
  83   TString parList ; //this will be list of parameters used for this analysis.\r
  84   const Int_t buffersize = 255;\r
  85   char onePar[buffersize] ;\r
  86   \r
  87   snprintf(onePar,buffersize,"--- AliAnaShowerParameter ---\n") ;\r
  88   parList+=onePar ;     \r
  89   snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;\r
  90   parList+=onePar ;\r
  91   snprintf(onePar,buffersize,"fNCellsCutMin: (Cut in the minimum number of cells) %e\n",fNCellsCutMin) ;\r
  92   parList+=onePar ;  \r
  93   snprintf(onePar,buffersize,"fNCellsCutMax: (Cut in the maximum number of cells) %e\n",fNCellsCutMax) ;\r
  94   parList+=onePar ;\r
  95   snprintf(onePar,buffersize,"fLambdaCut: (Cut in the minimum lambda) %e\n",fLambdaCut) ;\r
  96   parList+=onePar ;\r
  97   \r
  98   //Get parameters set in base class.\r
  99   parList += GetBaseParametersList() ;\r
 100   \r
 101   //Get parameters set in PID class.\r
 102   parList += GetCaloPID()->GetPIDParametersList() ;\r
 103   \r
 104   //Get parameters set in FiducialCut class (not available yet)\r
 105   //parlist += GetFidCut()->GetFidCutParametersList() \r
 106   \r
 107   return new TObjString(parList) ;\r
 108 }\r
 109 \r
 110 \r
 111 //________________________________________________________________________\r
 112 TList *  AliAnaShowerParameter::GetCreateOutputObjects()\r
 113 {  \r
 114   // Create histograms to be saved in output file and \r
 115   // store them in outputContainer\r
 116   TList * outputContainer = new TList() ; \r
 117   outputContainer->SetName("PhotonHistos") ; \r
 118         \r
 119   Int_t nptbins  = GetHistoPtBins();\r
 120   Int_t nphibins = GetHistoPhiBins();\r
 121   Int_t netabins = GetHistoEtaBins();\r
 122   Float_t ptmax  = GetHistoPtMax();\r
 123   Float_t phimax = GetHistoPhiMax();\r
 124   Float_t etamax = GetHistoEtaMax();\r
 125   Float_t ptmin  = GetHistoPtMin();\r
 126   Float_t phimin = GetHistoPhiMin();\r
 127   Float_t etamin = GetHistoEtaMin();    \r
 128   \r
 129   //General non-MC Cluster histograms\r
 130   fhNClusters = new TH1F ("hNClusters","NClusters",21,-0.5,20.5); \r
 131   fhNClusters->SetXTitle("N_{Clusters}");\r
 132   outputContainer->Add(fhNClusters) ;\r
 133   \r
 134   fhNCellCluster = new TH2F ("hNCellCluster","Number of cell per cluster",nptbins,ptmin,ptmax,21,-0.5,20.5); \r
 135   fhNCellCluster->SetYTitle("N_{Cells}");  \r
 136   fhNCellCluster->SetXTitle("p_{T}");\r
 137   outputContainer->Add(fhNCellCluster) ;\r
 138   \r
 139   fhEtaPhiPtCluster = new TH3F\r
 140   ("hEtaPhiPtCluster","#phi_{Cluster} and #eta_{Cluster}",nptbins,ptmin,ptmax,nphibins,phimin,phimax,netabins,etamin,etamax); \r
 141   fhEtaPhiPtCluster->SetZTitle("#eta");  \r
 142   fhEtaPhiPtCluster->SetYTitle("#phi");\r
 143   fhEtaPhiPtCluster->SetXTitle("pT_{Cluster} (GeV/c)");\r
 144   outputContainer->Add(fhEtaPhiPtCluster) ; \r
 145   \r
 146   fhLambdaPtCluster  = new TH2F\r
 147   ("hLambdaCluster","#lambda_{Cluster}",nptbins,ptmin,ptmax,300,0,3); \r
 148   fhLambdaPtCluster->SetYTitle("#lambda_{0}^{2}");\r
 149   fhLambdaPtCluster->SetXTitle("p_{T Cluster} (GeV/c)");\r
 150   outputContainer->Add(fhLambdaPtCluster) ;\r
 151   \r
 152   if(IsDataMC()){\r
 153     \r
 154     fhLambdaPtPhoton  = new TH2F\r
 155     ("hLambdaPtPhoton","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2); \r
 156     fhLambdaPtPhoton->SetYTitle("#lambda_{0}^{2}");\r
 157     fhLambdaPtPhoton->SetXTitle("pT_{#gamma, Reco} (GeV)");\r
 158     outputContainer->Add(fhLambdaPtPhoton) ;\r
 159     \r
 160     fhLambdaPtPi0  = new TH2F\r
 161     ("hLambdaPtPi0","#lambda_{#pi^{0}}",nptbins,ptmin,ptmax,200,0,2); \r
 162     fhLambdaPtPi0->SetYTitle("#lambda_{0}^{2}");\r
 163     fhLambdaPtPi0->SetXTitle("pT_{#pi^{0}, Reco} (GeV)");\r
 164     outputContainer->Add(fhLambdaPtPi0) ;\r
 165     \r
 166     fhLambdaPtPion  = new TH2F\r
 167     ("hLambdaPtPion","#lambda_{#pi^{+}}",nptbins,ptmin,ptmax,200,0,2); \r
 168     fhLambdaPtPion->SetYTitle("#lambda_{0}^{2}");\r
 169     fhLambdaPtPion->SetXTitle("pT_{#pi^{+}, Reco} (GeV)");\r
 170     outputContainer->Add(fhLambdaPtPion) ;\r
 171    \r
 172     fhPtTruthPi0  = new TH1D("hPtTruthPi0","#pi^{0} MC truth pT",nptbins,ptmin,ptmax) ;\r
 173     fhPtTruthPi0->SetXTitle("pT_{#pi^{0}, Reco} (GeV)");\r
 174     outputContainer->Add(fhPtTruthPi0) ;\r
 175     \r
 176   }//Histos with MC\r
 177   \r
 178   return outputContainer ;\r
 179   \r
 180 }\r
 181 \r
 182 //____________________________________________________________________________\r
 183 void AliAnaShowerParameter::Init()\r
 184 {\r
 185   \r
 186   //Init\r
 187   //Do some checks\r
 188   if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD()){\r
 189     printf("AliAnaShowerParameter::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");\r
 190     abort();\r
 191   }\r
 192   else  if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD()){\r
 193     printf("AliAnaShowerParameter::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");\r
 194     abort();\r
 195   }\r
 196   \r
 197 }\r
 198 \r
 199 //____________________________________________________________________________\r
 200 void AliAnaShowerParameter::InitParameters()\r
 201 {\r
 202   \r
 203   //Initialize the parameters of the analysis.\r
 204   AddToHistogramsName("AnaLambda_");\r
 205   \r
 206   fCalorimeter = "EMCAL" ;\r
 207         \r
 208   fTimeCutMin  = -1;\r
 209   fTimeCutMax  = 9999999;\r
 210   fNCellsCutMin = 0 ;\r
 211   fNCellsCutMax = 0 ;\r
 212   fLambdaCut = 0.01 ;\r
 213         \r
 214 }\r
 215 \r
 216 //__________________________________________________________________\r
 217 void  AliAnaShowerParameter::MakeAnalysisFillAOD() \r
 218 {\r
 219   /*//Do analysis and fill aods\r
 220   //Search for photons in fCalorimeter \r
 221   \r
 222   //Get vertex for photon momentum calculation\r
 223   \r
 224   for (Int_t iev = 0; iev < GetNMixedEvent(); iev++) {\r
 225     if (!GetMixedEvent()) \r
 226       GetReader()->GetVertex(GetVertex(iev));\r
 227     else \r
 228       GetMixedEvent()->GetVertexOfEvent(iev)->GetXYZ(GetVertex(iev)); \r
 229   } \r
 230   \r
 231   //Select the Calorimeter of the photon\r
 232   TObjArray * pl = 0x0; \r
 233   if(fCalorimeter == "PHOS")\r
 234     pl = GetPHOSClusters();\r
 235   else if (fCalorimeter == "EMCAL")\r
 236     pl = GetEMCALClusters();\r
 237   \r
 238   if(!pl){\r
 239     printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Careful cluster array NULL!!\n");\r
 240     return;\r
 241   }\r
 242   \r
 243   //Fill AODCaloClusters and AODParticle with PHOS/EMCAL aods\r
 244   TLorentzVector mom, mom2 ; \r
 245   Int_t nCaloClusters = pl->GetEntriesFast();   \r
 246   //Cut on the number of clusters in the event.\r
 247   if ((fNumClusters !=-1) && (nCaloClusters != fNumClusters)) return;\r
 248   Bool_t * indexConverted = new Bool_t[nCaloClusters];\r
 249   for (Int_t i = 0; i < nCaloClusters; i++) \r
 250     indexConverted[i] = kFALSE;\r
 251   \r
 252   for(Int_t icalo = 0; icalo < nCaloClusters; icalo++){    \r
 253     \r
 254     AliAODCaloCluster * calo =  (AliAODCaloCluster*) (pl->At(icalo));   \r
 255     Int_t evtIndex = 0 ; \r
 256     if (GetMixedEvent()) {\r
 257       evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ; \r
 258     }\r
 259     //Cluster selection, not charged, with photon id and in fiducial cut\r
 260           \r
 261     //Input from second AOD?\r
 262     Int_t input = 0;\r
 263     \r
 264     //Get Momentum vector, \r
 265     if (input == 0) \r
 266       calo->GetMomentum(mom,GetVertex(evtIndex)) ;//Assume that come from vertex in straight line\r
 267     \r
 268     //Skip the cluster if it doesn't fit inside the cuts.\r
 269     if(mom.Pt() < GetMinPt() || mom.Pt() > GetMaxPt() ) continue ; \r
 270     Double_t tof = calo->GetTOF()*1e9;    \r
 271     if(tof < fTimeCutMin || tof > fTimeCutMax) continue;          \r
 272     if(calo->GetNCells() <= fNCellsCut) continue;\r
 273     \r
 274     //Check acceptance selection\r
 275     if(IsFiducialCutOn()){\r
 276       Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;\r
 277       if(! in ) continue ;\r
 278     }\r
 279     \r
 280     //Create AOD for analysis\r
 281     AliAODPWG4Particle aodph = AliAODPWG4Particle(mom);\r
 282     Int_t label = calo->GetLabel();\r
 283     aodph.SetLabel(label);\r
 284     aodph.SetInputFileIndex(input);\r
 285     \r
 286     //Set the indices of the original caloclusters  \r
 287     aodph.SetCaloLabel(calo->GetID(),-1);\r
 288     aodph.SetDetector(fCalorimeter);\r
 289     if(GetDebug() > 1) \r
 290       printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Min pt cut and fiducial cut passed: pt %3.2f, phi %2.2f, eta %1.2f\n",aodph.Pt(),aodph.Phi(),aodph.Eta()); \r
 291     \r
 292     //Check Distance to Bad channel, set bit.\r
 293     Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel\r
 294     if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;\r
 295     if(distBad < fMinDist) //In bad channel (PHOS cristal size 2.2x2.2 cm)\r
 296       continue ;\r
 297     \r
 298     if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Bad channel cut passed %4.2f\n",distBad);\r
 299     \r
 300     if     (distBad > fMinDist3) aodph.SetDistToBad(2) ;\r
 301     else if(distBad > fMinDist2) aodph.SetDistToBad(1) ; \r
 302     else                         aodph.SetDistToBad(0) ;\r
 303     \r
 304     //Skip matched clusters with tracks\r
 305     if(fRejectTrackMatch && IsTrackMatched(calo)) continue ;\r
 306     if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - TrackMatching cut passed \n");\r
 307     \r
 308     //Set PID bits for later selection (AliAnaPi0 for example)\r
 309     //GetPDG already called in SetPIDBits.\r
 310     GetCaloPID()->SetPIDBits(fCalorimeter,calo,&aodph, GetCaloUtils());\r
 311     if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - PID Bits set \n");                    \r
 312     \r
 313     //Play with the MC stack if available\r
 314     //Check origin of the candidates\r
 315     if(IsDataMC()){\r
 316       aodph.SetTag(GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader(), aodph.GetInputFileIndex()));\r
 317       if(GetDebug() > 0) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Origin of candidate, bit map %d\n",aodph.GetTag());\r
 318     }\r
 319     \r
 320     // Check if cluster comes from a conversion in the material in front of the calorimeter\r
 321     // Do invariant mass of all pairs, if mass is close to 0, then it is conversion.\r
 322     \r
 323     if(fCheckConversion && nCaloClusters > 1){\r
 324       Bool_t bConverted = kFALSE;\r
 325       Int_t id2 = -1;\r
 326                   \r
 327       //Check if set previously as converted couple, if so skip its use.\r
 328       if (indexConverted[icalo]) continue;\r
 329                   \r
 330       for(Int_t jcalo = icalo + 1 ; jcalo < nCaloClusters ; jcalo++) {\r
 331         //Check if set previously as converted couple, if so skip its use.\r
 332         if (indexConverted[jcalo]) continue;\r
 333         //printf("Check Conversion indeces %d and %d\n",icalo,jcalo);\r
 334         AliAODCaloCluster * calo2 =  (AliAODCaloCluster*) (pl->At(jcalo));              //Get cluster kinematics\r
 335         Int_t evtIndex2 = 0 ; \r
 336         if (GetMixedEvent()) {\r
 337           evtIndex2=GetMixedEvent()->EventIndexForCaloCluster(calo2->GetID()) ; \r
 338         }        \r
 339         calo2->GetMomentum(mom2,GetVertex(evtIndex2));\r
 340         //Check only certain regions\r
 341         Bool_t in2 = kTRUE;\r
 342         if(IsFiducialCutOn()) in2 =  GetFiducialCut()->IsInFiducialCut(mom2,fCalorimeter) ;\r
 343         if(!in2) continue;      \r
 344         \r
 345         //Get mass of pair, if small, take this pair.\r
 346         //printf("\t both in calo, mass %f, cut %f\n",(mom+mom2).M(),fMassCut);\r
 347         if((mom+mom2).M() < fMassCut){  \r
 348           bConverted = kTRUE;\r
 349           id2 = calo2->GetID();\r
 350           indexConverted[jcalo]=kTRUE;\r
 351           break;\r
 352         }\r
 353                           \r
 354       }//Mass loop\r
 355                   \r
 356       if(bConverted){ \r
 357         if(fAddConvertedPairsToAOD){\r
 358           //Create AOD of pair analysis\r
 359           TLorentzVector mpair = mom+mom2;\r
 360           AliAODPWG4Particle aodpair = AliAODPWG4Particle(mpair);\r
 361           aodpair.SetLabel(aodph.GetLabel());\r
 362           aodpair.SetInputFileIndex(input);\r
 363           \r
 364           //printf("Index %d, Id %d\n",icalo, calo->GetID());\r
 365           //Set the indeces of the original caloclusters  \r
 366           aodpair.SetCaloLabel(calo->GetID(),id2);\r
 367           aodpair.SetDetector(fCalorimeter);\r
 368           aodpair.SetPdg(aodph.GetPdg());\r
 369           aodpair.SetTag(aodph.GetTag());\r
 370           \r
 371           //Add AOD with pair object to aod branch\r
 372           AddAODParticle(aodpair);\r
 373           //printf("\t \t both added pair\n");\r
 374         }\r
 375         \r
 376         //Do not add the current calocluster\r
 377         continue;\r
 378       }//converted pair\r
 379     }//check conversion\r
 380           \r
 381     //Add AOD with photon object to aod branch\r
 382     AddAODParticle(aodph);\r
 383     \r
 384   }//loop;\r
 385   delete [] indexConverted;\r
 386         \r
 387   if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD()  End fill AODs, with %d entries \n",GetOutputAODBranch()->GetEntriesFast());  \r
 388   \r
 389 */\r
 390 }\r
 391 //__________________________________________________________________\r
 392 void  AliAnaShowerParameter::MakeAnalysisFillHistograms() \r
 393 {\r
 394   \r
 395   //Do analysis and fill histograms\r
 396   if ( GetDebug() > 0) printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Starting analysis.\n");\r
 397     \r
 398   // Access MC information in stack if requested, check that it exists. \r
 399   AliStack * stack = 0x0;\r
 400   //TParticle * primary = 0x0;   \r
 401   TClonesArray * mcparticles0 = 0x0;\r
 402   //AliAODMCParticle * aodprimary = 0x0; \r
 403   //TObjArray * pl = 0x0;\r
 404   Int_t NClusters = 0 ;\r
 405     TLorentzVector momCluster ;\r
 406 \r
 407   \r
 408   //Check if the stack is available when analysing MC data.\r
 409   if(IsDataMC()){\r
 410     \r
 411     if(GetReader()->ReadStack()){\r
 412       stack =  GetMCStack() ;\r
 413       if(!stack) {\r
 414         printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Stack not available, is the MC handler called? STOP\n");\r
 415                                 abort();\r
 416       }\r
 417       \r
 418     }\r
 419     else if(GetReader()->ReadAODMCParticles()){\r
 420       \r
 421       //Get the list of MC particles\r
 422       mcparticles0 = GetReader()->GetAODMCParticles(0);\r
 423       if(!mcparticles0 && GetDebug() > 0)       {\r
 424         printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() -  Standard MCParticles not available !\n"); \r
 425       }\r
 426     }// is data and MC\r
 427   }     \r
 428   //Loop on stored AOD photons\r
 429 \r
 430   TClonesArray *  clustArray = 0x0 ; \r
 431   clustArray = GetAODCaloClusters() ;\r
 432   NClusters = clustArray->GetEntriesFast() ;    \r
 433   fhNClusters->Fill(NClusters) ;\r
 434 \r
 435   if(GetDebug() > 0) printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - aod branch entries %d\n", NClusters);\r
 436   \r
 437   for(Int_t iCluster = 0 ; iCluster < NClusters ; iCluster++){    \r
 438     Int_t input = 0 ;\r
 439 \r
 440     AliAODCaloCluster * caloCluster =  (AliAODCaloCluster*) (clustArray->At(iCluster)) ;\r
 441     caloCluster->GetMomentum(momCluster,GetVertex(0)) ;\r
 442     AliAODPWG4Particle * aodCluster = new AliAODPWG4Particle(momCluster) ;\r
 443     Int_t labelCluster = caloCluster->GetLabel() ;\r
 444     aodCluster->SetLabel(labelCluster) ;;\r
 445     aodCluster->SetInputFileIndex(input) ;\r
 446     aodCluster->SetCaloLabel(caloCluster->GetID(),-1) ;\r
 447     aodCluster->SetDetector(fCalorimeter) ;\r
 448     aodCluster->SetTag(GetMCAnalysisUtils()->CheckOrigin(caloCluster->GetLabels(),caloCluster->GetNLabels(),GetReader(), aodCluster->GetInputFileIndex())) ;\r
 449 \r
 450     if (GetDebug() > 3) printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Cluster particle label = %d\n",labelCluster) ;\r
 451     \r
 452     //Fill Cluster histograms \r
 453     Float_t ptcluster  = aodCluster->Pt() ;\r
 454     Float_t phicluster = aodCluster->Phi() ;\r
 455     Float_t etacluster = aodCluster->Eta() ;\r
 456     Float_t lambdaMainCluster = caloCluster->GetM02() ;\r
 457     Int_t iNumCell = caloCluster->GetNCells() ; \r
 458     \r
 459     if (iNumCell>=fNCellsCutMin&&iNumCell<=fNCellsCutMax&&lambdaMainCluster>=fLambdaCut){\r
 460       \r
 461       //Fill the basic non-MC information about the cluster.\r
 462       fhNCellCluster->Fill(ptcluster,iNumCell) ;\r
 463       fhEtaPhiPtCluster->Fill(ptcluster,phicluster,etacluster) ;\r
 464       fhLambdaPtCluster->Fill(ptcluster,lambdaMainCluster) ;\r
 465       \r
 466     //Play with the MC data if available\r
 467       if(IsDataMC()){\r
 468         \r
 469         if(GetReader()->ReadStack() && !stack) return;\r
 470         if(GetReader()->ReadAODMCParticles() && !mcparticles0) return;\r
 471         \r
 472         //Get the tag from AliMCAnalysisUtils for PID\r
 473         Int_t tag = aodCluster->GetTag();\r
 474         if (GetDebug() > 3) printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Cluster particle tag= %d\n",tag) ;\r
 475         if ( aodCluster->GetLabel() < 0){\r
 476           if(GetDebug() > 0) \r
 477             printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Label is -1; problem in the MC ESD? ");\r
 478           continue;\r
 479         }\r
 480         \r
 481         //Check if the tag matches one of the different particle types and fill the corresponding histograms\r
 482         if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)&&!(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0))) {\r
 483           fhLambdaPtPhoton->Fill(ptcluster,lambdaMainCluster) ;\r
 484         }//kMCPhoton\r
 485         \r
 486         if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) {     \r
 487           fhLambdaPtPi0->Fill(ptcluster,lambdaMainCluster) ;        \r
 488         }\r
 489         if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPion)) {\r
 490           fhLambdaPtPion->Fill(ptcluster,lambdaMainCluster) ;\r
 491         }\r
 492         \r
 493         // Access MC information in stack if requested, check that it exists.\r
 494         Int_t label = aodCluster->GetLabel();\r
 495         if(label < 0) {\r
 496           printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() *** bad label ***:  label %d \n", label);\r
 497           continue;\r
 498         }\r
 499       }\r
 500     }\r
 501   }\r
 502   \r
 503   Float_t fPtGener = 0 ;\r
 504   if(IsDataMC()){\r
 505     TParticle* pGener = stack->Particle(0) ;\r
 506     fPtGener = pGener->Pt() ;\r
 507   }\r
 508 \r
 509   \r
 510 }\r
 511 \r
 512 \r
 513 //__________________________________________________________________\r
 514 void AliAnaShowerParameter::Print(const Option_t * opt) const\r
 515 {\r
 516   //Print some relevant parameters set for the analysis\r
 517   \r
 518   if(! opt)\r
 519     return;\r
 520   \r
 521   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;\r
 522   AliAnaPartCorrBaseClass::Print(" ");\r
 523   printf("Calorimeter            =     %s\n", fCalorimeter.Data()) ;\r
 524   printf("Min number of cells in cluster is        > %f \n", fNCellsCutMin);\r
 525   printf("Max number of cells in cluster is        > %f \n", fNCellsCutMax);\r
 526   printf("Min lambda of cluster is        > %f \n", fLambdaCut);\r
 527   printf("    \n") ;\r
 528   \r
 529 } \r