PWGGA/CaloTrackCorrelations/AliAnaPhotonConvInCalo.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 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 //_________________________________________________________________________
  17 //
  18 // Conversions pairs analysis
  19 // Check if cluster comes from a conversion in the material in front of the calorimeter
  20 // Do invariant mass of all pairs, if mass is close to 0, then it is conversion.
  21 // Input are selected clusters with AliAnaPhoton
  22 //
  23 //
  24 //-- Author: Gustavo Conesa (LPSC-IN2P3-CNRS)
  25 //////////////////////////////////////////////////////////////////////////////
  26
  27
  28 // --- ROOT system ---
  29 #include <TH2F.h>
  30 #include <TH3D.h>
  31 #include <TClonesArray.h>
  32 #include <TObjString.h>
  33 #include "TParticle.h"
  34 #include "TDatabasePDG.h"
  35
  36 // --- Analysis system ---
  37 #include "AliAnaPhotonConvInCalo.h"
  38 #include "AliCaloTrackReader.h"
  39 #include "AliStack.h"
  40 #include "AliCaloPID.h"
  41 #include "AliMCAnalysisUtils.h"
  42 #include "AliFiducialCut.h"
  43 #include "AliVCluster.h"
  44 #include "AliAODMCParticle.h"
  45
  46 ClassImp(AliAnaPhotonConvInCalo)
  47
  48 //________________________________________________
  49 AliAnaPhotonConvInCalo::AliAnaPhotonConvInCalo() :
  50 AliAnaCaloTrackCorrBaseClass(),
  51 fRemoveConvertedPair(kFALSE),
  52 fAddConvertedPairsToAOD(kFALSE),
  53 fMassCut(0),
  54 fConvAsymCut(1.),                  fConvDEtaCut(2.),
  55 fConvDPhiMinCut(-1.),              fConvDPhiMaxCut(7.),
  56 fMomentum(),                       fProdVertex(),
  57 // Histograms
  58 fhPtPhotonConv(0),                 fhEtaPhiPhotonConv(0),          fhEtaPhi05PhotonConv(0),
  59 fhConvDeltaEta(0),                 fhConvDeltaPhi(0),              fhConvDeltaEtaPhi(0),
  60 fhConvAsym(0),                     fhConvPt(0),
  61 fhConvDistEta(0),                  fhConvDistEn(0),                fhConvDistMass(0),
  62 fhConvDistEtaCutEta(0),            fhConvDistEnCutEta(0),          fhConvDistMassCutEta(0),
  63 fhConvDistEtaCutMass(0),           fhConvDistEnCutMass(0),
  64 fhConvDistEtaCutAsy(0),            fhConvDistEnCutAsy(0),
  65
  66 // MC histograms
  67 fhPtConversionTagged(0),           fhPtAntiNeutronTagged(0),
  68 fhPtAntiProtonTagged(0),           fhPtUnknownTagged(0),
  69
  70 fhConvDeltaEtaMCConversion(0),     fhConvDeltaPhiMCConversion(0),  fhConvDeltaEtaPhiMCConversion(0),
  71 fhConvAsymMCConversion(0),         fhConvPtMCConversion(0),
  72 fhConvDispersionMCConversion(0),   fhConvM02MCConversion(0),
  73
  74 fhConvDeltaEtaMCAntiNeutron(0),    fhConvDeltaPhiMCAntiNeutron(0), fhConvDeltaEtaPhiMCAntiNeutron(0),
  75 fhConvAsymMCAntiNeutron(0),        fhConvPtMCAntiNeutron(0),
  76 fhConvDispersionMCAntiNeutron(0),  fhConvM02MCAntiNeutron(0),
  77 fhConvDeltaEtaMCAntiProton(0),     fhConvDeltaPhiMCAntiProton(0),  fhConvDeltaEtaPhiMCAntiProton(0),
  78 fhConvAsymMCAntiProton(0),         fhConvPtMCAntiProton(0),
  79 fhConvDispersionMCAntiProton(0),   fhConvM02MCAntiProton(0),
  80 fhConvDeltaEtaMCString(0),         fhConvDeltaPhiMCString(0),      fhConvDeltaEtaPhiMCString(0),
  81 fhConvAsymMCString(0),             fhConvPtMCString(0),
  82 fhConvDispersionMCString(0),       fhConvM02MCString(0),
  83 fhConvDistMCConversion(0),         fhConvDistMCConversionCuts(0)
  84 {
  85   //default ctor
  86
  87   //Initialize parameters
  88   InitParameters();
  89
  90 }
  91
  92 //_____________________________________________________
  93 TObjString *  AliAnaPhotonConvInCalo::GetAnalysisCuts()
  94 {
  95   //Save parameters used for analysis
  96   TString parList ; //this will be list of parameters used for this analysis.
  97   const Int_t buffersize = 255;
  98   char onePar[buffersize] ;
  99
 100   snprintf(onePar,buffersize,"--- AliAnaPhotonConvInCalo---:") ;
 101   parList+=onePar ;
 102   snprintf(onePar,buffersize,"Conversion Selection: fConvAsymCut %1.2f, fConvDEtaCut %1.2f fConvDPhiCut (%1.2f,%1.2f)",
 103            fConvAsymCut, fConvDEtaCut, fConvDPhiMinCut, fConvDPhiMaxCut) ;
 104   parList+=onePar ;
 105
 106   return new TObjString(parList) ;
 107 }
 108
 109 //_______________________________________________________
 110 TList *  AliAnaPhotonConvInCalo::GetCreateOutputObjects()
 111 {
 112   // Create histograms to be saved in output file and
 113   // store them in outputContainer
 114   TList * outputContainer = new TList() ;
 115   outputContainer->SetName("PhotonConvInCaloHistos") ;
 116
 117   Int_t nptbins  = GetHistogramRanges()->GetHistoPtBins();  Float_t ptmax  = GetHistogramRanges()->GetHistoPtMax();  Float_t ptmin  = GetHistogramRanges()->GetHistoPtMin();
 118   Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Float_t phimax = GetHistogramRanges()->GetHistoPhiMax(); Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();
 119   Int_t netabins = GetHistogramRanges()->GetHistoEtaBins(); Float_t etamax = GetHistogramRanges()->GetHistoEtaMax(); Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
 120
 121   fhPtPhotonConv  = new TH1F("hPtPhotonConv","Number of #gamma over calorimeter, conversion",nptbins,ptmin,ptmax);
 122   fhPtPhotonConv->SetYTitle("N");
 123   fhPtPhotonConv->SetXTitle("p_{T #gamma}(GeV/c)");
 124   outputContainer->Add(fhPtPhotonConv) ;
 125
 126   fhEtaPhiPhotonConv  = new TH2F
 127   ("hEtaPhiPhotonConv","#eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
 128   fhEtaPhiPhotonConv->SetYTitle("#phi (rad)");
 129   fhEtaPhiPhotonConv->SetXTitle("#eta");
 130   outputContainer->Add(fhEtaPhiPhotonConv) ;
 131   if(GetMinPt() < 0.5){
 132     fhEtaPhi05PhotonConv  = new TH2F
 133     ("hEtaPhi05PhotonConv","#eta vs #phi, E > 0.5",netabins,etamin,etamax,nphibins,phimin,phimax);
 134     fhEtaPhi05PhotonConv->SetYTitle("#phi (rad)");
 135     fhEtaPhi05PhotonConv->SetXTitle("#eta");
 136     outputContainer->Add(fhEtaPhi05PhotonConv) ;
 137   }
 138
 139   fhConvDeltaEta  = new TH2F
 140   ("hConvDeltaEta","#Delta #eta of selected conversion pairs",100,0,fMassCut,netabins*2,-0.5,0.5);
 141   fhConvDeltaEta->SetYTitle("#Delta #eta");
 142   fhConvDeltaEta->SetXTitle("Pair Mass (GeV/c^2)");
 143   outputContainer->Add(fhConvDeltaEta) ;
 144
 145   fhConvDeltaPhi  = new TH2F
 146   ("hConvDeltaPhi","#Delta #phi of selected conversion pairs",100,0,fMassCut,nphibins*2,-0.5,0.5);
 147   fhConvDeltaPhi->SetYTitle("#Delta #phi");
 148   fhConvDeltaPhi->SetXTitle("Pair Mass (GeV/c^2)");
 149   outputContainer->Add(fhConvDeltaPhi) ;
 150
 151   fhConvDeltaEtaPhi  = new TH2F
 152   ("hConvDeltaEtaPhi","#Delta #eta vs #Delta #phi of selected conversion pairs",netabins,-0.5,0.5,nphibins,-0.5,0.5);
 153   fhConvDeltaEtaPhi->SetYTitle("#Delta #phi");
 154   fhConvDeltaEtaPhi->SetXTitle("#Delta #eta");
 155   outputContainer->Add(fhConvDeltaEtaPhi) ;
 156
 157   fhConvAsym  = new TH2F
 158   ("hConvAsym","Asymmetry of selected conversion pairs",100,0,fMassCut,100,0,1);
 159   fhConvAsym->SetYTitle("Asymmetry");
 160   fhConvAsym->SetXTitle("Pair Mass (GeV/c^2)");
 161   outputContainer->Add(fhConvAsym) ;
 162
 163   fhConvPt  = new TH2F
 164   ("hConvPt","p_{T} of selected conversion pairs",100,0,fMassCut,100,0.,10.);
 165   fhConvPt->SetYTitle("Pair p_{T} (GeV/c)");
 166   fhConvPt->SetXTitle("Pair Mass (GeV/c^2)");
 167   outputContainer->Add(fhConvPt) ;
 168
 169   fhConvDistEta  = new TH2F
 170   ("hConvDistEta","distance to conversion vertex",100,-0.7,0.7,100,0.,5.);
 171   fhConvDistEta->SetXTitle("#eta");
 172   fhConvDistEta->SetYTitle(" distance (m)");
 173   outputContainer->Add(fhConvDistEta) ;
 174
 175   fhConvDistEn  = new TH2F
 176   ("hConvDistEn","distance to conversion vertex",nptbins,ptmin,ptmax,100,0.,5.);
 177   fhConvDistEn->SetXTitle("E (GeV)");
 178   fhConvDistEn->SetYTitle(" distance (m)");
 179   outputContainer->Add(fhConvDistEn) ;
 180
 181   fhConvDistMass  = new TH2F
 182   ("hConvDistMass","distance to conversion vertex",100,0,fMassCut,100,0.,5.);
 183   fhConvDistMass->SetXTitle("m (GeV/c^2)");
 184   fhConvDistMass->SetYTitle(" distance (m)");
 185   outputContainer->Add(fhConvDistMass) ;
 186
 187   fhConvDistEtaCutEta  = new TH2F
 188   ("hConvDistEtaCutEta","distance to conversion vertex, dEta < 0.05",100,-0.7,0.7,100,0.,5.);
 189   fhConvDistEtaCutEta->SetXTitle("#eta");
 190   fhConvDistEtaCutEta->SetYTitle(" distance (m)");
 191   outputContainer->Add(fhConvDistEtaCutEta) ;
 192
 193   fhConvDistEnCutEta  = new TH2F
 194   ("hConvDistEnCutEta","distance to conversion vertex, dEta < 0.05",nptbins,ptmin,ptmax,100,0.,5.);
 195   fhConvDistEnCutEta->SetXTitle("E (GeV)");
 196   fhConvDistEnCutEta->SetYTitle(" distance (m)");
 197   outputContainer->Add(fhConvDistEnCutEta) ;
 198
 199   fhConvDistMassCutEta  = new TH2F
 200   ("hConvDistMassCutEta","distance to conversion vertex, dEta < 0.05",100,0,fMassCut,100,0.,5.);
 201   fhConvDistMassCutEta->SetXTitle("m (GeV/c^2)");
 202   fhConvDistMassCutEta->SetYTitle(" distance (m)");
 203   outputContainer->Add(fhConvDistMassCutEta) ;
 204
 205   fhConvDistEtaCutMass  = new TH2F
 206   ("hConvDistEtaCutMass","distance to conversion vertex, dEta < 0.05, m < 10 MeV",100,-0.7,0.7,100,0.,5.);
 207   fhConvDistEtaCutMass->SetXTitle("#eta");
 208   fhConvDistEtaCutMass->SetYTitle(" distance (m)");
 209   outputContainer->Add(fhConvDistEtaCutMass) ;
 210
 211   fhConvDistEnCutMass  = new TH2F
 212   ("hConvDistEnCutMass","distance to conversion vertex, dEta < 0.05, m < 10 MeV",nptbins,ptmin,ptmax,100,0.,5.);
 213   fhConvDistEnCutMass->SetXTitle("E (GeV)");
 214   fhConvDistEnCutMass->SetYTitle(" distance (m)");
 215   outputContainer->Add(fhConvDistEnCutMass) ;
 216
 217   fhConvDistEtaCutAsy  = new TH2F
 218   ("hConvDistEtaCutAsy","distance to conversion vertex, dEta < 0.05, m < 10 MeV, A < 0.1",100,-0.7,0.7,100,0.,5.);
 219   fhConvDistEtaCutAsy->SetXTitle("#eta");
 220   fhConvDistEtaCutAsy->SetYTitle(" distance (m)");
 221   outputContainer->Add(fhConvDistEtaCutAsy) ;
 222
 223   fhConvDistEnCutAsy  = new TH2F
 224   ("hConvDistEnCutAsy","distance to conversion vertex, dEta < 0.05, m < 10 MeV, A < 0.1",nptbins,ptmin,ptmax,100,0.,5.);
 225   fhConvDistEnCutAsy->SetXTitle("E (GeV)");
 226   fhConvDistEnCutAsy->SetYTitle(" distance (m)");
 227   outputContainer->Add(fhConvDistEnCutAsy) ;
 228
 229   if(IsDataMC())
 230   {
 231     fhPtConversionTagged  = new TH1F("hPtMCConversionTagged","Number of converted #gamma over calorimeter, tagged as converted",nptbins,ptmin,ptmax);
 232     fhPtConversionTagged->SetYTitle("N");
 233     fhPtConversionTagged->SetXTitle("p_{T #gamma}(GeV/c)");
 234     outputContainer->Add(fhPtConversionTagged) ;
 235
 236
 237     fhPtAntiNeutronTagged  = new TH1F("hPtMCAntiNeutronTagged","Number of AntiNeutron id as Photon over calorimeter, tagged as converted",nptbins,ptmin,ptmax);
 238     fhPtAntiNeutronTagged->SetYTitle("N");
 239     fhPtAntiNeutronTagged->SetXTitle("p_{T #gamma}(GeV/c)");
 240     outputContainer->Add(fhPtAntiNeutronTagged) ;
 241
 242     fhPtAntiProtonTagged  = new TH1F("hPtMCAntiProtonTagged","Number of AntiProton id as Photon over calorimeter, tagged as converted",nptbins,ptmin,ptmax);
 243     fhPtAntiProtonTagged->SetYTitle("N");
 244     fhPtAntiProtonTagged->SetXTitle("p_{T #gamma}(GeV/c)");
 245     outputContainer->Add(fhPtAntiProtonTagged) ;
 246
 247     fhPtUnknownTagged  = new TH1F("hPtMCUnknownTagged","Number of Unknown id as Photon over calorimeter, tagged as converted",nptbins,ptmin,ptmax);
 248     fhPtUnknownTagged->SetYTitle("N");
 249     fhPtUnknownTagged->SetXTitle("p_{T #gamma}(GeV/c)");
 250     outputContainer->Add(fhPtUnknownTagged) ;
 251
 252     fhConvDeltaEtaMCConversion  = new TH2F
 253     ("hConvDeltaEtaMCConversion","#Delta #eta of selected conversion pairs from real conversions",100,0,fMassCut,netabins,-0.5,0.5);
 254     fhConvDeltaEtaMCConversion->SetYTitle("#Delta #eta");
 255     fhConvDeltaEtaMCConversion->SetXTitle("Pair Mass (GeV/c^2)");
 256     outputContainer->Add(fhConvDeltaEtaMCConversion) ;
 257
 258     fhConvDeltaPhiMCConversion  = new TH2F
 259     ("hConvDeltaPhiMCConversion","#Delta #phi of selected conversion pairs from real conversions",100,0,fMassCut,nphibins,-0.5,0.5);
 260     fhConvDeltaPhiMCConversion->SetYTitle("#Delta #phi");
 261     fhConvDeltaPhiMCConversion->SetXTitle("Pair Mass (GeV/c^2)");
 262     outputContainer->Add(fhConvDeltaPhiMCConversion) ;
 263
 264     fhConvDeltaEtaPhiMCConversion  = new TH2F
 265     ("hConvDeltaEtaPhiMCConversion","#Delta #eta vs #Delta #phi of selected conversion pairs, from real conversions",netabins,-0.5,0.5,nphibins,-0.5,0.5);
 266     fhConvDeltaEtaPhiMCConversion->SetYTitle("#Delta #phi");
 267     fhConvDeltaEtaPhiMCConversion->SetXTitle("#Delta #eta");
 268     outputContainer->Add(fhConvDeltaEtaPhiMCConversion) ;
 269
 270     fhConvAsymMCConversion  = new TH2F
 271     ("hConvAsymMCConversion","Asymmetry of selected conversion pairs from real conversions",100,0,fMassCut,100,0,1);
 272     fhConvAsymMCConversion->SetYTitle("Asymmetry");
 273     fhConvAsymMCConversion->SetXTitle("Pair Mass (GeV/c^2)");
 274     outputContainer->Add(fhConvAsymMCConversion) ;
 275
 276     fhConvPtMCConversion  = new TH2F
 277     ("hConvPtMCConversion","p_{T} of selected conversion pairs from real conversions",100,0,fMassCut,100,0.,10.);
 278     fhConvPtMCConversion->SetYTitle("Pair p_{T} (GeV/c)");
 279     fhConvPtMCConversion->SetXTitle("Pair Mass (GeV/c^2)");
 280     outputContainer->Add(fhConvPtMCConversion) ;
 281
 282     fhConvDispersionMCConversion  = new TH2F
 283     ("hConvDispersionMCConversion","p_{T} of selected conversion pairs from real conversions",100,0.,1.,100,0.,1.);
 284     fhConvDispersionMCConversion->SetYTitle("Dispersion cluster 1");
 285     fhConvDispersionMCConversion->SetXTitle("Dispersion cluster 2");
 286     outputContainer->Add(fhConvDispersionMCConversion) ;
 287
 288     fhConvM02MCConversion  = new TH2F
 289     ("hConvM02MCConversion","p_{T} of selected conversion pairs from string",100,0.,1.,100,0.,1.);
 290     fhConvM02MCConversion->SetYTitle("M02 cluster 1");
 291     fhConvM02MCConversion->SetXTitle("M02 cluster 2");
 292     outputContainer->Add(fhConvM02MCConversion) ;
 293
 294     fhConvDeltaEtaMCAntiNeutron  = new TH2F
 295     ("hConvDeltaEtaMCAntiNeutron","#Delta #eta of selected conversion pairs from anti-neutrons",100,0,fMassCut,netabins,-0.5,0.5);
 296     fhConvDeltaEtaMCAntiNeutron->SetYTitle("#Delta #eta");
 297     fhConvDeltaEtaMCAntiNeutron->SetXTitle("Pair Mass (GeV/c^2)");
 298     outputContainer->Add(fhConvDeltaEtaMCAntiNeutron) ;
 299
 300     fhConvDeltaPhiMCAntiNeutron  = new TH2F
 301     ("hConvDeltaPhiMCAntiNeutron","#Delta #phi of selected conversion pairs from anti-neutrons",100,0,fMassCut,nphibins,-0.5,0.5);
 302     fhConvDeltaPhiMCAntiNeutron->SetYTitle("#Delta #phi");
 303     fhConvDeltaPhiMCAntiNeutron->SetXTitle("Pair Mass (GeV/c^2)");
 304     outputContainer->Add(fhConvDeltaPhiMCAntiNeutron) ;
 305
 306     fhConvDeltaEtaPhiMCAntiNeutron  = new TH2F
 307     ("hConvDeltaEtaPhiMCAntiNeutron","#Delta #eta vs #Delta #phi of selected conversion pairs from anti-neutrons",netabins,-0.5,0.5,nphibins,-0.5,0.5);
 308     fhConvDeltaEtaPhiMCAntiNeutron->SetYTitle("#Delta #phi");
 309     fhConvDeltaEtaPhiMCAntiNeutron->SetXTitle("#Delta #eta");
 310     outputContainer->Add(fhConvDeltaEtaPhiMCAntiNeutron) ;
 311
 312     fhConvAsymMCAntiNeutron  = new TH2F
 313     ("hConvAsymMCAntiNeutron","Asymmetry of selected conversion pairs from anti-neutrons",100,0,fMassCut,100,0,1);
 314     fhConvAsymMCAntiNeutron->SetYTitle("Asymmetry");
 315     fhConvAsymMCAntiNeutron->SetXTitle("Pair Mass (GeV/c^2)");
 316     outputContainer->Add(fhConvAsymMCAntiNeutron) ;
 317
 318     fhConvPtMCAntiNeutron  = new TH2F
 319     ("hConvPtMCAntiNeutron","p_{T} of selected conversion pairs from anti-neutrons",100,0,fMassCut,100,0.,10.);
 320     fhConvPtMCAntiNeutron->SetYTitle("Pair p_{T} (GeV/c)");
 321     fhConvPtMCAntiNeutron->SetXTitle("Pair Mass (GeV/c^2)");
 322     outputContainer->Add(fhConvPtMCAntiNeutron) ;
 323
 324     fhConvDispersionMCAntiNeutron  = new TH2F
 325     ("hConvDispersionMCAntiNeutron","p_{T} of selected conversion pairs from anti-neutrons",100,0.,1.,100,0.,1.);
 326     fhConvDispersionMCAntiNeutron->SetYTitle("Dispersion cluster 1");
 327     fhConvDispersionMCAntiNeutron->SetXTitle("Dispersion cluster 2");
 328     outputContainer->Add(fhConvDispersionMCAntiNeutron) ;
 329
 330     fhConvM02MCAntiNeutron  = new TH2F
 331     ("hConvM02MCAntiNeutron","p_{T} of selected conversion pairs from string",100,0.,1.,100,0.,1.);
 332     fhConvM02MCAntiNeutron->SetYTitle("M02 cluster 1");
 333     fhConvM02MCAntiNeutron->SetXTitle("M02 cluster 2");
 334     outputContainer->Add(fhConvM02MCAntiNeutron) ;
 335
 336     fhConvDeltaEtaMCAntiProton  = new TH2F
 337     ("hConvDeltaEtaMCAntiProton","#Delta #eta of selected conversion pairs from anti-protons",100,0,fMassCut,netabins,-0.5,0.5);
 338     fhConvDeltaEtaMCAntiProton->SetYTitle("#Delta #eta");
 339     fhConvDeltaEtaMCAntiProton->SetXTitle("Pair Mass (GeV/c^2)");
 340     outputContainer->Add(fhConvDeltaEtaMCAntiProton) ;
 341
 342     fhConvDeltaPhiMCAntiProton  = new TH2F
 343     ("hConvDeltaPhiMCAntiProton","#Delta #phi of selected conversion pairs from anti-protons",100,0,fMassCut,nphibins,-0.5,0.5);
 344     fhConvDeltaPhiMCAntiProton->SetYTitle("#Delta #phi");
 345     fhConvDeltaPhiMCAntiProton->SetXTitle("Pair Mass (GeV/c^2)");
 346     outputContainer->Add(fhConvDeltaPhiMCAntiProton) ;
 347
 348     fhConvDeltaEtaPhiMCAntiProton  = new TH2F
 349     ("hConvDeltaEtaPhiMCAntiProton","#Delta #eta vs #Delta #phi of selected conversion pairs from anti-protons",netabins,-0.5,0.5,nphibins,-0.5,0.5);
 350     fhConvDeltaEtaPhiMCAntiProton->SetYTitle("#Delta #phi");
 351     fhConvDeltaEtaPhiMCAntiProton->SetXTitle("#Delta #eta");
 352     outputContainer->Add(fhConvDeltaEtaPhiMCAntiProton) ;
 353
 354     fhConvAsymMCAntiProton  = new TH2F
 355     ("hConvAsymMCAntiProton","Asymmetry of selected conversion pairs from anti-protons",100,0,fMassCut,100,0,1);
 356     fhConvAsymMCAntiProton->SetYTitle("Asymmetry");
 357     fhConvAsymMCAntiProton->SetXTitle("Pair Mass (GeV/c^2)");
 358     outputContainer->Add(fhConvAsymMCAntiProton) ;
 359
 360     fhConvPtMCAntiProton  = new TH2F
 361     ("hConvPtMCAntiProton","p_{T} of selected conversion pairs from anti-protons",100,0,fMassCut,100,0.,10.);
 362     fhConvPtMCAntiProton->SetYTitle("Pair p_{T} (GeV/c)");
 363     fhConvPtMCAntiProton->SetXTitle("Pair Mass (GeV/c^2)");
 364     outputContainer->Add(fhConvPtMCAntiProton) ;
 365
 366     fhConvDispersionMCAntiProton  = new TH2F
 367     ("hConvDispersionMCAntiProton","p_{T} of selected conversion pairs from anti-protons",100,0.,1.,100,0.,1.);
 368     fhConvDispersionMCAntiProton->SetYTitle("Dispersion cluster 1");
 369     fhConvDispersionMCAntiProton->SetXTitle("Dispersion cluster 2");
 370     outputContainer->Add(fhConvDispersionMCAntiProton) ;
 371
 372     fhConvM02MCAntiProton  = new TH2F
 373     ("hConvM02MCAntiProton","p_{T} of selected conversion pairs from string",100,0.,1.,100,0.,1.);
 374     fhConvM02MCAntiProton->SetYTitle("M02 cluster 1");
 375     fhConvM02MCAntiProton->SetXTitle("M02 cluster 2");
 376     outputContainer->Add(fhConvM02MCAntiProton) ;
 377
 378     fhConvDeltaEtaMCString  = new TH2F
 379     ("hConvDeltaEtaMCString","#Delta #eta of selected conversion pairs from string",100,0,fMassCut,netabins,-0.5,0.5);
 380     fhConvDeltaEtaMCString->SetYTitle("#Delta #eta");
 381     fhConvDeltaEtaMCString->SetXTitle("Pair Mass (GeV/c^2)");
 382     outputContainer->Add(fhConvDeltaEtaMCString) ;
 383
 384     fhConvDeltaPhiMCString  = new TH2F
 385     ("hConvDeltaPhiMCString","#Delta #phi of selected conversion pairs from string",100,0,fMassCut,nphibins,-0.5,0.5);
 386     fhConvDeltaPhiMCString->SetYTitle("#Delta #phi");
 387     fhConvDeltaPhiMCString->SetXTitle("Pair Mass (GeV/c^2)");
 388     outputContainer->Add(fhConvDeltaPhiMCString) ;
 389
 390     fhConvDeltaEtaPhiMCString  = new TH2F
 391     ("hConvDeltaEtaPhiMCString","#Delta #eta vs #Delta #phi of selected conversion pairs from string",netabins,-0.5,0.5,nphibins,-0.5,0.5);
 392     fhConvDeltaEtaPhiMCString->SetYTitle("#Delta #phi");
 393     fhConvDeltaEtaPhiMCString->SetXTitle("#Delta #eta");
 394     outputContainer->Add(fhConvDeltaEtaPhiMCString) ;
 395
 396     fhConvAsymMCString  = new TH2F
 397     ("hConvAsymMCString","Asymmetry of selected conversion pairs from string",100,0,fMassCut,100,0,1);
 398     fhConvAsymMCString->SetYTitle("Asymmetry");
 399     fhConvAsymMCString->SetXTitle("Pair Mass (GeV/c^2)");
 400     outputContainer->Add(fhConvAsymMCString) ;
 401
 402     fhConvPtMCString  = new TH2F
 403     ("hConvPtMCString","p_{T} of selected conversion pairs from string",100,0,fMassCut,100,0.,10.);
 404     fhConvPtMCString->SetYTitle("Pair p_{T} (GeV/c)");
 405     fhConvPtMCString->SetXTitle("Pair Mass (GeV/c^2)");
 406     outputContainer->Add(fhConvPtMCString) ;
 407
 408     fhConvDispersionMCString  = new TH2F
 409     ("hConvDispersionMCString","p_{T} of selected conversion pairs from string",100,0.,1.,100,0.,1.);
 410     fhConvDispersionMCString->SetYTitle("Dispersion cluster 1");
 411     fhConvDispersionMCString->SetXTitle("Dispersion cluster 2");
 412     outputContainer->Add(fhConvDispersionMCString) ;
 413
 414     fhConvM02MCString  = new TH2F
 415     ("hConvM02MCString","p_{T} of selected conversion pairs from string",100,0.,1.,100,0.,1.);
 416     fhConvM02MCString->SetYTitle("M02 cluster 1");
 417     fhConvM02MCString->SetXTitle("M02 cluster 2");
 418     outputContainer->Add(fhConvM02MCString) ;
 419
 420     fhConvDistMCConversion  = new TH2F
 421     ("hConvDistMCConversion","calculated conversion distance vs real vertes for MC conversion",100,0.,5.,100,0.,5.);
 422     fhConvDistMCConversion->SetYTitle("distance");
 423     fhConvDistMCConversion->SetXTitle("vertex R");
 424     outputContainer->Add(fhConvDistMCConversion) ;
 425
 426     fhConvDistMCConversionCuts  = new TH2F
 427     ("hConvDistMCConversionCuts","calculated conversion distance vs real vertes for MC conversion, deta < 0.05, m < 10 MeV, asym < 0.1",100,0.,5.,100,0.,5.);
 428     fhConvDistMCConversionCuts->SetYTitle("distance");
 429     fhConvDistMCConversionCuts->SetXTitle("vertex R");
 430     outputContainer->Add(fhConvDistMCConversionCuts) ;
 431
 432   }
 433
 434   return outputContainer ;
 435
 436 }
 437
 438 //___________________________________________
 439 void AliAnaPhotonConvInCalo::InitParameters()
 440 {
 441   //Initialize the parameters of the analysis.
 442   AddToHistogramsName("AnaPhotonConvInCalo_");
 443
 444   fMassCut                = 0.03; //30 MeV
 445   fRemoveConvertedPair    = kFALSE;
 446   fAddConvertedPairsToAOD = kFALSE;
 447
 448 }
 449
 450 //_________________________________________________
 451 void  AliAnaPhotonConvInCalo::MakeAnalysisFillAOD()
 452 {
 453   //Do conversion photon analysis and fill aods
 454
 455   //Loop on stored AOD photons
 456   Int_t naod = GetOutputAODBranch()->GetEntriesFast();
 457   AliDebug(1,Form("AOD branch entries %d", naod));
 458
 459   //List to be used in conversion analysis, to tag the cluster as candidate for conversion
 460   Bool_t * indexConverted = new Bool_t[naod];
 461   for (Int_t i = 0; i < naod; i++) indexConverted[i] = kFALSE;
 462
 463   for(Int_t iaod = 0; iaod < naod ; iaod++)
 464   {
 465     AliAODPWG4Particle* calo =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
 466
 467     Bool_t bConverted = kFALSE;
 468     Int_t id2 = -1;
 469
 470     //Check if set previously as converted couple, if so skip its use.
 471     if (indexConverted[iaod]) continue;
 472
 473     // Second cluster loop
 474     AliAODPWG4Particle* calo2 = 0;
 475     for(Int_t jaod = iaod + 1 ; jaod < naod ; jaod++)
 476     {
 477       //Check if set previously as converted couple, if so skip its use.
 478       if (indexConverted[jaod]) continue;
 479       //printf("Check Conversion indeces %d and %d\n",iaod,jaod);
 480       calo2 =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(jaod));
 481
 482       //................................................
 483       //Get mass of pair, if small, take this pair.
 484       Float_t pairM     = calo->GetPairMass(calo2);
 485       //printf("\t both in calo, mass %f, cut %f\n",pairM,fMassCut);
 486       if(pairM < fMassCut)
 487       {
 488         calo->SetTagged(kFALSE);
 489         id2 = calo2->GetCaloLabel(0);
 490         Float_t asymmetry = TMath::Abs(calo->E()-calo2->E())/(calo->E()+calo2->E());
 491         Float_t dPhi      = (calo->Momentum())->Phi()-(calo2->Momentum())->Phi();
 492         Float_t dEta      = (calo->Momentum())->Eta()-(calo2->Momentum())->Eta();
 493
 494         //...............................................
 495         //Fill few histograms with kinematics of the pair
 496         //FIXME, move all this to MakeAnalysisFillHistograms ...
 497
 498         fhConvDeltaEta   ->Fill( pairM, dPhi      );
 499         fhConvDeltaPhi   ->Fill( pairM, dEta      );
 500         fhConvAsym       ->Fill( pairM, asymmetry );
 501         fhConvDeltaEtaPhi->Fill( dEta , dPhi      );
 502         fhConvPt         ->Fill( pairM, (calo->Momentum())->Pt()+(calo2->Momentum())->Pt());
 503
 504         //Estimate conversion distance, T. Awes, M. Ivanov
 505         //Under the assumption that the pair has zero mass, and that each electron
 506         //of the pair has the same momentum, they will each have the same bend radius
 507         //given by R=p/(qB) = p / (300 B) with p in [MeV/c], B in [Tesla] and R in [m].
 508         //With nominal ALICE magnet current of 30kA B=0.5T, and so with E_cluster=p,
 509         //R = E/1.5 [cm].  Under these assumptions, the distance from the conversion
 510         //point to the MCEal can be related to the separation distance, L=2y, on the MCEal
 511         //as d = sqrt(R^2 -(R-y)^2) = sqrt(2Ry - y^2). And since R>>y we can write as
 512         //d = sqrt(E*L/1.5) where E is the cluster energy and L is the distance in cm between
 513         //the clusters.
 514
 515         TObjArray * clusters    = 0;
 516         if(calo->GetDetectorTag() == kEMCAL) clusters = GetEMCALClusters();
 517         else                                 clusters = GetPHOSClusters();
 518
 519
 520         Int_t iclus = -1;
 521         AliVCluster *cluster1 = FindCluster(clusters,calo ->GetCaloLabel(0),iclus);
 522         AliVCluster *cluster2 = FindCluster(clusters,calo2->GetCaloLabel(0),iclus);
 523
 524         Float_t pos1[3];
 525         cluster1->GetPosition(pos1);
 526         Float_t pos2[3];
 527         cluster2->GetPosition(pos2);
 528         Float_t clustDist = TMath::Sqrt((pos1[0]-pos2[0])*(pos1[0]-pos2[0])+
 529                                         (pos1[1]-pos2[1])*(pos1[1]-pos2[1])+
 530                                         (pos1[2]-pos2[2])*(pos1[2]-pos2[2]));
 531
 532         Float_t convDist  = TMath::Sqrt(calo->E() *clustDist*0.01/0.15);
 533         Float_t convDist2 = TMath::Sqrt(calo2->E()*clustDist*0.01/0.15);
 534         //printf("l = %f, e1 = %f, d1=%f, e2 = %f, d2=%f\n",clustDist,calo->E(),convDist,calo2->E(),convDist2);
 535         AliDebug(2,Form("Pair with mass %2.3f < %2.3f, %1.2f < dPhi %2.2f < %2.2f, dEta %f < %2.2f, asymmetry %2.2f< %2.2f; \n    cluster1 id %d, e %2.3f  SM %d, eta %2.3f, phi %2.3f ; \n    cluster2 id %d, e %2.3f, SM %d,eta %2.3f, phi %2.3f",
 536                         pairM,fMassCut,fConvDPhiMinCut, dPhi, fConvDPhiMaxCut, dEta, fConvDEtaCut, asymmetry, fConvAsymCut,
 537                         calo->GetCaloLabel(0),calo->E(),GetCaloUtils()->GetModuleNumber(calo,GetReader()->GetInputEvent()), calo->Eta(), calo->Phi(),
 538                         id2, calo2->E(), GetCaloUtils()->GetModuleNumber(calo2,GetReader()->GetInputEvent()),calo2->Eta(), calo2->Phi()));
 539
 540         fhConvDistEta ->Fill(calo ->Eta(),convDist );
 541         fhConvDistEta ->Fill(calo2->Eta(),convDist2);
 542         fhConvDistEn  ->Fill(calo ->E(),  convDist );
 543         fhConvDistEn  ->Fill(calo2->E(),  convDist2);
 544         fhConvDistMass->Fill(pairM, convDist );
 545         //dEta cut
 546         if(dEta<0.05)
 547         {
 548           fhConvDistEtaCutEta ->Fill(calo->Eta(), convDist );
 549           fhConvDistEtaCutEta ->Fill(calo2->Eta(),convDist2);
 550           fhConvDistEnCutEta  ->Fill(calo->E(),   convDist );
 551           fhConvDistEnCutEta  ->Fill(calo2->E(),  convDist2);
 552           fhConvDistMassCutEta->Fill(pairM, convDist );
 553           //mass cut
 554           if(pairM<0.01)
 555           {//10 MeV
 556             fhConvDistEtaCutMass ->Fill(calo ->Eta(), convDist );
 557             fhConvDistEtaCutMass ->Fill(calo2->Eta(), convDist2);
 558             fhConvDistEnCutMass  ->Fill(calo ->E(),   convDist );
 559             fhConvDistEnCutMass  ->Fill(calo2->E(),   convDist2);
 560             // asymmetry cut
 561             if(asymmetry<0.1)
 562             {
 563               fhConvDistEtaCutAsy ->Fill(calo ->Eta(), convDist );
 564               fhConvDistEtaCutAsy ->Fill(calo2->Eta(), convDist2);
 565               fhConvDistEnCutAsy  ->Fill(calo ->E(),   convDist );
 566               fhConvDistEnCutAsy  ->Fill(calo2->E(),   convDist2);
 567             }//asymmetry cut
 568           }//mass cut
 569         }//dEta cut
 570
 571         //...............................................
 572         //Select pairs in a eta-phi window
 573         if(TMath::Abs(dEta) < fConvDEtaCut    &&
 574            TMath::Abs(dPhi) < fConvDPhiMaxCut &&
 575            TMath::Abs(dPhi) > fConvDPhiMinCut &&
 576            asymmetry        < fConvAsymCut       )
 577         {
 578           indexConverted[iaod] = kTRUE;
 579           indexConverted[jaod] = kTRUE;
 580           bConverted           = kTRUE;
 581         }
 582         //printf("Accepted? %d\n",bConverted);
 583         //...........................................
 584         //Fill more histograms, simulated data
 585         //FIXME, move all this to MakeAnalysisFillHistograms ...
 586         if(IsDataMC())
 587         {
 588           //Check the origin of the pair, look for conversion, antinucleons or jet correlations (strings)
 589           Int_t ancPDG    = 0;
 590           Int_t ancStatus = 0;
 591           Int_t ancLabel  = GetMCAnalysisUtils()->CheckCommonAncestor(cluster1->GetLabel(), cluster2->GetLabel(),
 592                                                                       GetReader(), ancPDG, ancStatus, fMomentum, fProdVertex);
 593
 594           // printf("AliAnaPhotonConvInCalo::MakeAnalysisFillHistograms() - Common ancestor label %d, pdg %d, name %s, status %d; \n",
 595           //                          ancLabel,ancPDG,TDatabasePDG::Instance()->GetParticle(ancPDG)->GetName(),ancStatus);
 596
 597           Int_t tag1 = calo ->GetTag();
 598           Int_t tag2 = calo2->GetTag();
 599           if(GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCConversion))
 600           {
 601             if(GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCConversion) && (ancPDG==22 || TMath::Abs(ancPDG)==11) && ancLabel > -1)
 602             {
 603               fhConvDeltaEtaMCConversion   ->Fill( pairM, dEta      );
 604               fhConvDeltaPhiMCConversion   ->Fill( pairM, dPhi      );
 605               fhConvAsymMCConversion       ->Fill( pairM, asymmetry );
 606               fhConvDeltaEtaPhiMCConversion->Fill( dEta , dPhi      );
 607               fhConvPtMCConversion         ->Fill( pairM, calo->Pt()+calo2->Pt());
 608               fhConvDispersionMCConversion ->Fill( cluster1->GetDispersion(), cluster2->GetDispersion());
 609               fhConvM02MCConversion        ->Fill( cluster1->GetM02(), cluster2->GetM02());
 610               fhConvDistMCConversion       ->Fill( convDist , fProdVertex.Mag() );
 611               fhConvDistMCConversion       ->Fill( convDist2, fProdVertex.Mag() );
 612
 613               if(dEta<0.05 && pairM<0.01 && asymmetry<0.1)
 614               {
 615                 fhConvDistMCConversionCuts->Fill( convDist , fProdVertex.Mag() );
 616                 fhConvDistMCConversionCuts->Fill( convDist2, fProdVertex.Mag() );
 617               }
 618
 619             }
 620           }
 621           else if(GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCAntiNeutron))
 622           {
 623             if(GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCAntiNeutron) && ancPDG==-2112 && ancLabel > -1)
 624             {
 625               fhConvDeltaEtaMCAntiNeutron    ->Fill( pairM, dEta      );
 626               fhConvDeltaPhiMCAntiNeutron    ->Fill( pairM, dPhi      );
 627               fhConvAsymMCAntiNeutron        ->Fill( pairM, asymmetry );
 628               fhConvDeltaEtaPhiMCAntiNeutron ->Fill( dEta , dPhi      );
 629               fhConvPtMCAntiNeutron          ->Fill( pairM, calo->Pt()+calo2->Pt());
 630               fhConvDispersionMCAntiNeutron  ->Fill( cluster1->GetDispersion(), cluster2->GetDispersion());
 631               fhConvM02MCAntiNeutron         ->Fill( cluster1->GetM02(), cluster2->GetM02());
 632             }
 633           }
 634           else if(GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCAntiProton))
 635           {
 636             if(GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCAntiProton) && ancPDG==-2212 && ancLabel > -1)
 637             {
 638               fhConvDeltaEtaMCAntiProton    ->Fill( pairM, dEta      );
 639               fhConvDeltaPhiMCAntiProton    ->Fill( pairM, dPhi      );
 640               fhConvAsymMCAntiProton        ->Fill( pairM, asymmetry );
 641               fhConvDeltaEtaPhiMCAntiProton ->Fill( dEta , dPhi      );
 642               fhConvPtMCAntiProton          ->Fill( pairM, calo->Pt()+calo2->Pt());
 643               fhConvDispersionMCAntiProton  ->Fill( cluster1->GetDispersion(), cluster2->GetDispersion());
 644               fhConvM02MCAntiProton         ->Fill( cluster1->GetM02(), cluster2->GetM02());
 645             }
 646           }
 647
 648           //Pairs coming from fragmenting pairs.
 649           if( ancPDG < 22 && ancLabel > 7 && (ancStatus == 11 || ancStatus == 12) )
 650           {
 651             fhConvDeltaEtaMCString    ->Fill( pairM, dPhi);
 652             fhConvDeltaPhiMCString    ->Fill( pairM, dPhi);
 653             fhConvAsymMCString        ->Fill( pairM, TMath::Abs(calo->E()-calo2->E())/(calo->E()+calo2->E()) );
 654             fhConvDeltaEtaPhiMCString ->Fill( dPhi,  dPhi);
 655             fhConvPtMCString          ->Fill( pairM, calo->Pt()+calo2->Pt());
 656             fhConvDispersionMCString  ->Fill( cluster1->GetDispersion(), cluster2->GetDispersion());
 657             fhConvM02MCString         ->Fill( cluster1->GetM02(), cluster2->GetM02());
 658           }
 659
 660         }// Data MC
 661
 662         break;
 663       }
 664
 665     }//Mass loop
 666
 667     //..........................................................................................................
 668     //Pair selected as converted, remove both clusters or recombine them into a photon and put them in the AOD
 669     if(bConverted)
 670     {
 671       //Add to AOD
 672       if(fAddConvertedPairsToAOD)
 673       {
 674         //Create AOD of pair analysis
 675         fMomentum = *(calo->Momentum())+*(calo2->Momentum());
 676         AliAODPWG4Particle aodpair = AliAODPWG4Particle(fMomentum);
 677         aodpair.SetLabel(calo->GetLabel());
 678
 679         //printf("Index %d, Id %d\n",iaod, calo->GetID());
 680         //Set the indeces of the original caloclusters
 681         aodpair.SetCaloLabel(calo->GetCaloLabel(0),id2);
 682         aodpair.SetDetectorTag(calo->GetDetectorTag());
 683         aodpair.SetIdentifiedParticleType(calo->GetIdentifiedParticleType());
 684         aodpair.SetTag(calo ->GetTag());
 685         aodpair.SetTagged(kTRUE);
 686         //Add AOD with pair object to aod branch
 687         AddAODParticle(aodpair);
 688         //printf("\t \t both added pair\n");
 689       }
 690
 691       //Do not add the current calocluster
 692       if(!fRemoveConvertedPair)
 693       {
 694         //printf("TAGGED\n");
 695         //Tag this cluster as likely conversion
 696         calo->SetTagged(kTRUE);
 697       }
 698     }//converted pair
 699
 700   }// main loop
 701
 702   // Remove entries identified as conversion electrons
 703   // Revise if this is OK
 704   if(fRemoveConvertedPair || fAddConvertedPairsToAOD)
 705   {
 706     for(Int_t iaod = 0; iaod < naod ; iaod++)
 707     {
 708       if(indexConverted[iaod])GetOutputAODBranch()->RemoveAt(iaod);
 709     }
 710     GetOutputAODBranch()->Compress();
 711   }
 712
 713   delete [] indexConverted;
 714
 715   AliDebug(1,Form("End fill AODs, with %d entries",GetOutputAODBranch()->GetEntriesFast()));
 716
 717 }
 718
 719 //________________________________________________________
 720 void  AliAnaPhotonConvInCalo::MakeAnalysisFillHistograms()
 721 {
 722   //Fill histograms
 723
 724   //-------------------------------------------------------------------
 725   // Access MC information in stack if requested, check that it exists.
 726   AliStack         * stack       = 0x0;
 727   TParticle        * primary     = 0x0;
 728   TClonesArray     * mcparticles = 0x0;
 729   AliAODMCParticle * aodprimary  = 0x0;
 730
 731   if(IsDataMC())
 732   {
 733     if(GetReader()->ReadStack())
 734     {
 735       stack =  GetMCStack() ;
 736       if(!stack)
 737       {
 738         AliFatal("Stack not available, is the MC handler called? STOP");
 739         return;
 740       }
 741     }
 742     else if(GetReader()->ReadAODMCParticles())
 743     {
 744       //Get the list of MC particles
 745       mcparticles = GetReader()->GetAODMCParticles();
 746       if(!mcparticles)
 747       {
 748         AliFatal("Standard MCParticles not available!");
 749         return;
 750       }
 751     }
 752   }// is data and MC
 753
 754   //----------------------------------
 755   //Loop on stored AOD photons
 756   Int_t naod = GetOutputAODBranch()->GetEntriesFast();
 757   AliDebug(1,Form("AOD branch entries %d", naod));
 758
 759   for(Int_t iaod = 0; iaod < naod ; iaod++)
 760   {
 761     AliAODPWG4Particle* ph =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
 762     //Int_t pdg = ph->GetIdentifiedParticleType();
 763
 764     if(ph->IsTagged())
 765     {
 766       AliDebug(2,Form("ID Photon: pt %f, phi %f, eta %f", ph->Pt(),ph->Phi(),ph->Eta())) ;
 767       //................................
 768       //Fill photon histograms
 769       Float_t ptcluster  = ph->Pt();
 770       Float_t phicluster = ph->Phi();
 771       Float_t etacluster = ph->Eta();
 772       Float_t ecluster   = ph->E();
 773
 774       fhPtPhotonConv->Fill(ptcluster);
 775       if(ecluster > 0.5)        fhEtaPhiPhotonConv  ->Fill(etacluster, phicluster);
 776       else if(GetMinPt() < 0.5) fhEtaPhi05PhotonConv->Fill(etacluster, phicluster);
 777
 778       //.......................................
 779       //Play with the MC data if available
 780       if(IsDataMC())
 781       {
 782         //....................................................................
 783         // Access MC information in stack if requested, check that it exists.
 784         Int_t label =ph->GetLabel();
 785         if(label < 0)
 786         {
 787           AliDebug(1,Form("*** bad label ***:  label %d", label));
 788           continue;
 789         }
 790
 791         //Float_t eprim   = 0;
 792         //Float_t ptprim  = 0;
 793         if(GetReader()->ReadStack())
 794         {
 795           if(label >=  stack->GetNtrack())
 796           {
 797             AliDebug(1,Form("*** large label ***:  label %d, n tracks %d", label, stack->GetNtrack()));
 798             continue ;
 799           }
 800
 801           primary = stack->Particle(label);
 802           if(!primary)
 803           {
 804             AliDebug(1,Form("*** no primary ***:  label %d", label));
 805             continue;
 806           }
 807           //eprim   = primary->Energy();
 808           //ptprim  = primary->Pt();
 809
 810         }
 811         else if(GetReader()->ReadAODMCParticles())
 812         {
 813           if(label >=  mcparticles->GetEntriesFast())
 814           {
 815             AliDebug(2,Form("*** large label ***:  label %d, n tracks %d",label, mcparticles->GetEntriesFast()));
 816             continue ;
 817           }
 818
 819           //Get the particle
 820           aodprimary = (AliAODMCParticle*) mcparticles->At(label);
 821
 822           if(!aodprimary)
 823           {
 824             AliDebug(2,Form("*** no primary ***:  label %d", label));
 825             continue;
 826           }
 827
 828           //eprim   = aodprimary->E();
 829           //ptprim  = aodprimary->Pt();
 830
 831         }
 832
 833         Int_t tag =ph->GetTag();
 834
 835         if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
 836         {
 837           if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion))
 838           {
 839             fhPtConversionTagged ->Fill(ptcluster);
 840           }
 841         }
 842         else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron))
 843         {
 844           fhPtAntiNeutronTagged ->Fill(ptcluster);
 845         }
 846         else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton))
 847         {
 848           fhPtAntiProtonTagged ->Fill(ptcluster);
 849         }
 850         else
 851         {
 852           fhPtUnknownTagged ->Fill(ptcluster);
 853         }
 854
 855       }//Histograms with MC
 856     }// tagged by conversion
 857   }// aod loop
 858
 859 }
 860
 861
 862 //________________________________________________________
 863 void AliAnaPhotonConvInCalo::Print(const Option_t * opt) const
 864 {
 865   //Print some relevant parameters set for the analysis
 866
 867   if(! opt)
 868     return;
 869
 870   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
 871   AliAnaCaloTrackCorrBaseClass::Print(" ");
 872
 873   printf("Add conversion pair to AOD           = %d\n",fAddConvertedPairsToAOD);
 874   printf("Conversion pair mass cut             = %f\n",fMassCut);
 875   printf("Conversion selection cut : A < %1.2f; %1.3f < Dphi < %1.3f; Deta < %1.3f\n",
 876          fConvAsymCut,fConvDPhiMinCut, fConvDPhiMaxCut, fConvDEtaCut);
 877
 878   printf("    \n") ;
 879
 880 }