PWGGA/PHOSTasks/PHOS_embedding/AliAnalysisTaskPi0DiffEfficiency.cxx

   1 #include "TChain.h"
   2 #include "TTree.h"
   3 #include "TObjArray.h"
   4 #include "TF1.h"
   5 #include "TFile.h"
   6 #include "TH1F.h"
   7 #include "TH2F.h"
   8 #include "TH2I.h"
   9 #include "TH3F.h"
  10 #include "TParticle.h"
  11 #include "TCanvas.h"
  12 #include "TStyle.h"
  13 #include "TRandom.h"
  14
  15 #include "AliAODMCParticle.h"
  16 #include "AliAnalysisManager.h"
  17 #include "AliMCEventHandler.h"
  18 #include "AliMCEvent.h"
  19 #include "AliStack.h"
  20 #include "AliAnalysisTaskSE.h"
  21 #include "AliAnalysisTaskPi0DiffEfficiency.h"
  22 #include "AliCaloPhoton.h"
  23 #include "AliPHOSGeometry.h"
  24 #include "AliPHOSAodCluster.h"
  25 #include "AliPHOSCalibData.h"
  26 #include "AliAODEvent.h"
  27 #include "AliAODCaloCluster.h"
  28 #include "AliAODVertex.h"
  29 #include "AliESDtrackCuts.h"
  30 #include "AliLog.h"
  31 #include "AliPID.h"
  32 #include "AliCDBManager.h"
  33 #include "AliCentrality.h"
  34
  35 // Analysis task to calculate pi0 registration efficiency
  36 // as a difference between spectra before and after embedding
  37 // Authors: Dmitri Peressounko
  38 // Date   : Aug.2011
  39
  40 ClassImp(AliAnalysisTaskPi0DiffEfficiency)
  41
  42 //________________________________________________________________________
  43 AliAnalysisTaskPi0DiffEfficiency::AliAnalysisTaskPi0DiffEfficiency(const char *name)
  44 : AliAnalysisTaskPi0Efficiency(name),
  45   fPHOSEvent1(0),
  46   fPHOSEvent2(0)
  47 {
  48
  49
  50
  51 }
  52
  53 //________________________________________________________________________
  54 void AliAnalysisTaskPi0DiffEfficiency::UserCreateOutputObjects()
  55 {
  56   // Create histograms
  57   // Called once
  58
  59   // ESD histograms
  60   if(fOutputContainer != NULL){
  61     delete fOutputContainer;
  62   }
  63   fOutputContainer = new TList();
  64   fOutputContainer->SetOwner(kTRUE);
  65
  66   //Event selection
  67   fOutputContainer->Add(new TH1F("hSelEvents","Event celection", 10,0.,10.)) ;
  68
  69   //vertex distribution
  70   fOutputContainer->Add(new TH1F("hZvertex","Z vertex position", 50,-25.,25.)) ;
  71
  72   //Centrality
  73   fOutputContainer->Add(new TH1F("hCentrality","Event centrality", 100,0.,100.)) ;
  74
  75   //QA histograms
  76   fOutputContainer->Add(new TH2F("hCluM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
  77   fOutputContainer->Add(new TH2F("hCluM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
  78   fOutputContainer->Add(new TH2F("hCluM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
  79
  80   Int_t nM       = 500;
  81   Double_t mMin  = 0.0;
  82   Double_t mMax  = 1.0;
  83   Int_t nPt      = 200;
  84   Double_t ptMin = 0;
  85   Double_t ptMax = 20;
  86
  87   char key[55] ;
  88   for(Int_t cent=0; cent<6; cent++){
  89     //Single photon
  90
  91     fOutputContainer->Add(new TH1F(Form("hPhotAll_DistBad2_cen%d",cent),"All clusters",nPt,ptMin,ptMax));
  92     fOutputContainer->Add(new TH1F(Form("hPhotAll_DistBad4_cen%d",cent),"All clusters",nPt,ptMin,ptMax));
  93     fOutputContainer->Add(new TH1F(Form("hPhotAll_DistBad6_cen%d",cent),"All clusters",nPt,ptMin,ptMax));
  94     snprintf(key,55,"hPhotAll_cen%d",cent) ;
  95     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
  96     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
  97     snprintf(key,55,"hPhotAllcore_cen%d",cent) ;
  98     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
  99     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 100     snprintf(key,55,"hPhotAllwou_cen%d",cent) ;
 101     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 102     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 103     snprintf(key,55,"hPhotCPV_cen%d",cent) ;
 104     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 105     snprintf(key,55,"hPhotCPVcore_cen%d",cent) ;
 106     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 107     snprintf(key,55,"hPhotCPV2_cen%d",cent) ;
 108     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 109     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 110     snprintf(key,55,"hPhotDisp_cen%d",cent) ;
 111     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 112     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 113     snprintf(key,55,"hPhotDisp2_cen%d",cent) ;
 114     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 115     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 116     snprintf(key,55,"hPhotDispcore_cen%d",cent) ;
 117     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 118     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 119     snprintf(key,55,"hPhotBoth_cen%d",cent) ;
 120     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 121     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 122     snprintf(key,55,"hPhotBothcore_cen%d",cent) ;
 123     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 124     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 125
 126     snprintf(key,55,"hNegPhotAll_cen%d",cent) ;
 127     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 128     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 129     snprintf(key,55,"hNegPhotAllcore_cen%d",cent) ;
 130     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 131     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 132     snprintf(key,55,"hNegPhotAllwou_cen%d",cent) ;
 133     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 134     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 135     snprintf(key,55,"hNegPhotCPV_cen%d",cent) ;
 136     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 137     snprintf(key,55,"hNegPhotCPVcore_cen%d",cent) ;
 138     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 139     snprintf(key,55,"hNegPhotCPV2_cen%d",cent) ;
 140     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 141     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 142     snprintf(key,55,"hNegPhotDisp_cen%d",cent) ;
 143     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 144     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 145     snprintf(key,55,"hNegPhotDisp2_cen%d",cent) ;
 146     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 147     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 148     snprintf(key,55,"hNegPhotDispcore_cen%d",cent) ;
 149     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 150     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 151     snprintf(key,55,"hNegPhotBoth_cen%d",cent) ;
 152     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 153     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 154     snprintf(key,55,"hNegPhotBothcore_cen%d",cent) ;
 155     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 156     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 157
 158     snprintf(key,55,"hOldMassPtAll_cen%d",cent) ;
 159     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 160     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 161     snprintf(key,55,"hOldMassPtAllcore_cen%d",cent) ;
 162     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 163     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 164     snprintf(key,55,"hOldMassPtAllwou_cen%d",cent) ;
 165     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 166     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 167     snprintf(key,55,"hOldMassPtCPV_cen%d",cent) ;
 168     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 169     snprintf(key,55,"hOldMassPtCPVcore_cen%d",cent) ;
 170     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 171     snprintf(key,55,"hOldMassPtCPV2_cen%d",cent) ;
 172     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 173     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 174     snprintf(key,55,"hOldMassPtDisp_cen%d",cent) ;
 175     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 176     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 177     snprintf(key,55,"hOldMassPtDisp2_cen%d",cent) ;
 178     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 179     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 180     snprintf(key,55,"hOldMassPtDispcore_cen%d",cent) ;
 181     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 182     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 183     snprintf(key,55,"hOldMassPtBoth_cen%d",cent) ;
 184     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 185     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 186     snprintf(key,55,"hOldMassPtBothcore_cen%d",cent) ;
 187     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 188     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 189
 190     snprintf(key,55,"hNewMassPtAll_cen%d",cent) ;
 191     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 192     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 193     snprintf(key,55,"hNewMassPtAllcore_cen%d",cent) ;
 194     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 195     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 196     snprintf(key,55,"hNewMassPtAllwou_cen%d",cent) ;
 197     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 198     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 199     snprintf(key,55,"hNewMassPtCPV_cen%d",cent) ;
 200     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 201     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 202     snprintf(key,55,"hNewMassPtCPVcore_cen%d",cent) ;
 203     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 204     snprintf(key,55,"hNewMassPtCPV2_cen%d",cent) ;
 205     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 206     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 207     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 208     snprintf(key,55,"hNewMassPtDisp_cen%d",cent) ;
 209     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 210     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 211     snprintf(key,55,"hNewMassPtDisp2_cen%d",cent) ;
 212     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 213     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 214     snprintf(key,55,"hNewMassPtDispcore_cen%d",cent) ;
 215     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 216     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 217     snprintf(key,55,"hNewMassPtBoth_cen%d",cent) ;
 218     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 219     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 220     snprintf(key,55,"hNewMassPtBothcore_cen%d",cent) ;
 221     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 222     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 223
 224     snprintf(key,55,"hNegMassPtAll_cen%d",cent) ;
 225     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 226     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 227     snprintf(key,55,"hNegMassPtAllcore_cen%d",cent) ;
 228     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 229     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 230     snprintf(key,55,"hNegMassPtAllwou_cen%d",cent) ;
 231     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 232     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 233     snprintf(key,55,"hNegMassPtCPV_cen%d",cent) ;
 234     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 235     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 236     snprintf(key,55,"hNegMassPtCPVcore_cen%d",cent) ;
 237     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 238     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 239     snprintf(key,55,"hNegMassPtCPV2_cen%d",cent) ;
 240     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 241     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 242     snprintf(key,55,"hNegMassPtDisp_cen%d",cent) ;
 243     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 244     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 245     snprintf(key,55,"hNegMassPtDisp2_cen%d",cent) ;
 246     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 247     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 248     snprintf(key,55,"hNegMassPtDispcore_cen%d",cent) ;
 249     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 250     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 251     snprintf(key,55,"hNegMassPtBoth_cen%d",cent) ;
 252     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 253     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 254     snprintf(key,55,"hNegMassPtBothcore_cen%d",cent) ;
 255     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 256     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 257
 258
 259     snprintf(key,55,"hMassPtAll_cen%d",cent) ;
 260     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 261     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 262     snprintf(key,55,"hMassPtAllcore_cen%d",cent) ;
 263     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 264     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 265     snprintf(key,55,"hMassPtAllwou_cen%d",cent) ;
 266     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 267     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 268     snprintf(key,55,"hMassPtCPV_cen%d",cent) ;
 269     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 270     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 271     snprintf(key,55,"hMassPtCPVcore_cen%d",cent) ;
 272     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 273     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 274     snprintf(key,55,"hMassPtCPV2_cen%d",cent) ;
 275     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 276     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 277     snprintf(key,55,"hMassPtDisp_cen%d",cent) ;
 278     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 279     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 280     snprintf(key,55,"hMassPtDisp2_cen%d",cent) ;
 281     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 282     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 283     snprintf(key,55,"hMassPtDispcore_cen%d",cent) ;
 284     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 285     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 286     snprintf(key,55,"hMassPtBoth_cen%d",cent) ;
 287     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 288     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 289     snprintf(key,55,"hMassPtBothcore_cen%d",cent) ;
 290     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 291     ((TH2F*)fOutputContainer->Last())->Sumw2() ;
 292
 293     snprintf(key,55,"hMassPtAll_a07_cen%d",cent) ;
 294     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 295     snprintf(key,55,"hMassPtCPV_a07_cen%d",cent) ;
 296     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 297     snprintf(key,55,"hMassPtCPV2_a07_cen%d",cent) ;
 298     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 299     snprintf(key,55,"hMassPtDisp_a07_cen%d",cent) ;
 300     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 301     snprintf(key,55,"hMassPtBoth_a07_cen%d",cent) ;
 302     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 303
 304     snprintf(key,55,"hMassPtAll_a08_cen%d",cent) ;
 305     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 306     snprintf(key,55,"hMassPtCPV_a08_cen%d",cent) ;
 307     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 308     snprintf(key,55,"hMassPtCPV2_a08_cen%d",cent) ;
 309     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 310     snprintf(key,55,"hMassPtDisp_a08_cen%d",cent) ;
 311     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 312     snprintf(key,55,"hMassPtBoth_a08_cen%d",cent) ;
 313     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 314
 315     snprintf(key,55,"hMassPtAll_a09_cen%d",cent) ;
 316     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 317     snprintf(key,55,"hMassPtCPV_a09_cen%d",cent) ;
 318     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 319     snprintf(key,55,"hMassPtCPV2_a09_cen%d",cent) ;
 320     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 321     snprintf(key,55,"hMassPtDisp_a09_cen%d",cent) ;
 322     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 323     snprintf(key,55,"hMassPtBoth_a09_cen%d",cent) ;
 324     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 325
 326
 327     //Mixed
 328     snprintf(key,55,"hMiMassPtAll_cen%d",cent) ;
 329     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 330     snprintf(key,55,"hMiMassPtAllcore_cen%d",cent) ;
 331     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 332     snprintf(key,55,"hMiMassPtAllwou_cen%d",cent) ;
 333     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 334     snprintf(key,55,"hMiMassPtCPV_cen%d",cent) ;
 335     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 336     snprintf(key,55,"hMiMassPtCPVcore_cen%d",cent) ;
 337     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 338     snprintf(key,55,"hMiMassPtCPV2_cen%d",cent) ;
 339     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 340     snprintf(key,55,"hMiMassPtDisp_cen%d",cent) ;
 341     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 342     snprintf(key,55,"hMiMassPtDisp2_cen%d",cent) ;
 343     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 344     snprintf(key,55,"hMiMassPtDispcore_cen%d",cent) ;
 345     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 346     snprintf(key,55,"hMiMassPtBoth_cen%d",cent) ;
 347     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 348     snprintf(key,55,"hMiMassPtBothcore_cen%d",cent) ;
 349     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 350
 351     snprintf(key,55,"hMiMassPtAll_a07_cen%d",cent) ;
 352     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 353     snprintf(key,55,"hMiMassPtCPV_a07_cen%d",cent) ;
 354     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 355     snprintf(key,55,"hMiMassPtCPV2_a07_cen%d",cent) ;
 356     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 357     snprintf(key,55,"hMiMassPtDisp_a07_cen%d",cent) ;
 358     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 359     snprintf(key,55,"hMiMassPtBoth_a07_cen%d",cent) ;
 360     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 361
 362     snprintf(key,55,"hMiMassPtAll_a08_cen%d",cent) ;
 363     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 364     snprintf(key,55,"hMiMassPtCPV_a08_cen%d",cent) ;
 365     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 366     snprintf(key,55,"hMiMassPtCPV2_a08_cen%d",cent) ;
 367     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 368     snprintf(key,55,"hMiMassPtDisp_a08_cen%d",cent) ;
 369     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 370     snprintf(key,55,"hMiMassPtBoth_a08_cen%d",cent) ;
 371     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 372
 373      snprintf(key,55,"hMiMassPtAll_a09_cen%d",cent) ;
 374     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 375     snprintf(key,55,"hMiMassPtCPV_a09_cen%d",cent) ;
 376     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 377     snprintf(key,55,"hMiMassPtCPV2_a09_cen%d",cent) ;
 378     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 379     snprintf(key,55,"hMiMassPtDisp_a09_cen%d",cent) ;
 380     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 381     snprintf(key,55,"hMiMassPtBoth_a09_cen%d",cent) ;
 382     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 383
 384     snprintf(key,55,"hMCMassPtAll_cen%d",cent) ;
 385     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 386     snprintf(key,55,"hMCMassPtAllcore_cen%d",cent) ;
 387     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 388     snprintf(key,55,"hMCMassPtAllwou_cen%d",cent) ;
 389     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 390     snprintf(key,55,"hMCMassPtCPV_cen%d",cent) ;
 391     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 392     snprintf(key,55,"hMCMassPtCPVcore_cen%d",cent) ;
 393     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 394     snprintf(key,55,"hMCMassPtCPV2_cen%d",cent) ;
 395     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 396     snprintf(key,55,"hMCMassPtDisp_cen%d",cent) ;
 397     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 398     snprintf(key,55,"hMCMassPtDisp2_cen%d",cent) ;
 399     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 400     snprintf(key,55,"hMCMassPtDispcore_cen%d",cent) ;
 401     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 402     snprintf(key,55,"hMCMassPtBoth_cen%d",cent) ;
 403     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 404     snprintf(key,55,"hMCMassPtBothcore_cen%d",cent) ;
 405     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 406
 407     snprintf(key,55,"hMCMassPtAll_a07_cen%d",cent) ;
 408     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 409     snprintf(key,55,"hMCMassPtCPV_a07_cen%d",cent) ;
 410     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 411     snprintf(key,55,"hMCMassPtCPV2_a07_cen%d",cent) ;
 412     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 413     snprintf(key,55,"hMCMassPtDisp_a07_cen%d",cent) ;
 414     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 415     snprintf(key,55,"hMCMassPtBoth_a07_cen%d",cent) ;
 416     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 417
 418     snprintf(key,55,"hMCMassPtAll_a08_cen%d",cent) ;
 419     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 420     snprintf(key,55,"hMCMassPtCPV_a08_cen%d",cent) ;
 421     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 422     snprintf(key,55,"hMCMassPtCPV2_a08_cen%d",cent) ;
 423     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 424     snprintf(key,55,"hMCMassPtDisp_a08_cen%d",cent) ;
 425     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 426     snprintf(key,55,"hMCMassPtBoth_a08_cen%d",cent) ;
 427     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 428
 429     snprintf(key,55,"hMCMassPtAll_a09_cen%d",cent) ;
 430     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 431     snprintf(key,55,"hMCMassPtCPV_a09_cen%d",cent) ;
 432     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 433     snprintf(key,55,"hMCMassPtCPV2_a09_cen%d",cent) ;
 434     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 435     snprintf(key,55,"hMCMassPtDisp_a09_cen%d",cent) ;
 436     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 437     snprintf(key,55,"hMCMassPtBoth_a09_cen%d",cent) ;
 438     fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 439
 440     //Single photon
 441     snprintf(key,55,"hMCPhotAll_cen%d",cent) ;
 442     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 443     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 444     snprintf(key,55,"hMCPhotAllcore_cen%d",cent) ;
 445     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 446     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 447     snprintf(key,55,"hMCPhotAllwou_cen%d",cent) ;
 448     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 449     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 450     snprintf(key,55,"hMCPhotCPV_cen%d",cent) ;
 451     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 452     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 453     snprintf(key,55,"hMCPhotCPVcore_cen%d",cent) ;
 454     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 455     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 456     snprintf(key,55,"hMCPhotCPV2_cen%d",cent) ;
 457     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 458     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 459     snprintf(key,55,"hMCPhotDisp_cen%d",cent) ;
 460     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 461     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 462     snprintf(key,55,"hMCPhotDisp2_cen%d",cent) ;
 463     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 464     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 465     snprintf(key,55,"hMCPhotDispcore_cen%d",cent) ;
 466     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 467     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 468     snprintf(key,55,"hMCPhotBoth_cen%d",cent) ;
 469     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 470     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 471     snprintf(key,55,"hMCPhotBothcore_cen%d",cent) ;
 472     fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
 473     ((TH1F*)fOutputContainer->Last())->Sumw2() ;
 474
 475   }
 476
 477   fOutputContainer->Add(new TH2F("hMCPi0M11","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 478   fOutputContainer->Add(new TH2F("hMCPi0M22","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 479   fOutputContainer->Add(new TH2F("hMCPi0M33","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 480   fOutputContainer->Add(new TH2F("hMCPi0M12","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 481   fOutputContainer->Add(new TH2F("hMCPi0M13","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 482   fOutputContainer->Add(new TH2F("hMCPi0M23","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
 483
 484   //MC
 485   for(Int_t cent=0; cent<6; cent++){
 486     snprintf(key,55,"hMC_rap_gamma_cen%d",cent) ;
 487     fOutputContainer->Add(new TH1F(key,"Rapidity pi0",200,-1.,1.)) ;
 488     snprintf(key,55,"hMC_rap_pi0_cen%d",cent) ;
 489     fOutputContainer->Add(new TH1F(key,"Rapidity pi0",200,-1.,1.)) ;
 490     snprintf(key,55,"hMC_rap_eta_cen%d",cent) ;
 491     fOutputContainer->Add(new TH1F("hMC_rap_eta","Rapidity eta",200,-1.,1.)) ;
 492     snprintf(key,55,"hMC_phi_gamma_cen%d",cent) ;
 493     fOutputContainer->Add(new TH1F(key,"Phi pi0",200,0.,TMath::TwoPi())) ;
 494     snprintf(key,55,"hMC_phi_pi0_cen%d",cent) ;
 495     fOutputContainer->Add(new TH1F(key,"Phi pi0",200,0.,TMath::TwoPi())) ;
 496     snprintf(key,55,"hMC_phi_eta_cen%d",cent) ;
 497     fOutputContainer->Add(new TH1F(key,"Phi eta",200,0.,TMath::TwoPi())) ;
 498     snprintf(key,55,"hMC_all_gamma_cen%d",cent) ;
 499     fOutputContainer->Add(new TH1F(key,"Rapidity photon",250,0.,25.)) ;
 500     snprintf(key,55,"hMC_all_pi0_cen%d",cent) ;
 501     fOutputContainer->Add(new TH1F(key,"Rapidity pi0",250,0.,25.)) ;
 502     snprintf(key,55,"hMC_all_eta_cen%d",cent) ;
 503     fOutputContainer->Add(new TH1F(key,"Rapidity eta",250,0.,25.)) ;
 504     snprintf(key,55,"hMC_unitEta_gamma_cen%d",cent) ;
 505     fOutputContainer->Add(new TH1F(key,"Pt photon",250,0.,25.)) ;
 506     snprintf(key,55,"hMC_unitEta_pi0_cen%d",cent) ;
 507     fOutputContainer->Add(new TH1F(key,"Rapidity eta",250,0.,25.)) ;
 508     snprintf(key,55,"hMC_unitEta_eta_cen%d",cent) ;
 509     fOutputContainer->Add(new TH1F(key,"Rapidity eta",250,0.,25.)) ;
 510   }
 511
 512   PostData(1, fOutputContainer);
 513
 514 }
 515
 516 //________________________________________________________________________
 517 void AliAnalysisTaskPi0DiffEfficiency::UserExec(Option_t *)
 518 {
 519   // Main loop, called for each event
 520   // Analyze ESD/AOD
 521
 522   FillHistogram("hSelEvents",0.5) ;
 523
 524   AliAODEvent *event = dynamic_cast<AliAODEvent*>(InputEvent());
 525   if (!event) {
 526     Printf("ERROR: Could not retrieve event");
 527     PostData(1, fOutputContainer);
 528     return;
 529   }
 530
 531   FillHistogram("hSelEvents",1.5) ;
 532   AliAODHeader *header = event->GetHeader() ;
 533
 534   // Checks if we have a primary vertex
 535   // Get primary vertices form ESD
 536   const AliAODVertex *esdVertex5 = event->GetPrimaryVertex();
 537
 538  // don't rely on ESD vertex, assume (0,0,0)
 539   Double_t vtx0[3] ={0.,0.,0.};
 540
 541
 542   FillHistogram("hZvertex",esdVertex5->GetZ());
 543   if (TMath::Abs(esdVertex5->GetZ()) > 10. ){
 544     PostData(1, fOutputContainer);
 545     return;
 546   }
 547   FillHistogram("hSelEvents",2.5) ;
 548
 549   //Vtx class z-bin
 550   //  Int_t zvtx = (Int_t)((vtx5[2]+10.)/2.) ;
 551   //  if(zvtx<0)zvtx=0 ;
 552   //  if(zvtx>9)zvtx=9 ;
 553   Int_t zvtx=0 ;
 554
 555 //  fCentrality=header->GetCentralityP()->GetCentralityPercentile("V0M"); // returns the centrality percentile,
 556 //                                                          //a float from 0 to 100 (or to the trigger efficiency)
 557    fCentrality=header->GetZDCN2Energy() ;
 558
 559   if( fCentrality < 0. || fCentrality>80.){
 560     PostData(1, fOutputContainer);
 561     return;
 562   }
 563   FillHistogram("hSelEvents",3.5) ;
 564   Float_t bins[7]={0.,5.,10.,20.,40.,60.,80.} ;
 565   fCenBin=0 ;
 566   while(fCenBin<6 && fCentrality > bins[fCenBin+1])
 567     fCenBin++ ;
 568
 569
 570   //reaction plain
 571   fRPfull= header->GetZDCN1Energy() ;
 572   if(fRPfull==999){ //reaction plain was not defined
 573     PostData(1, fOutputContainer);
 574     return;
 575   }
 576
 577   FillHistogram("hSelEvents",4.5) ;
 578   //All event selections done
 579   FillHistogram("hCentrality",fCentrality) ;
 580   //Reaction plain is defined in the range (-pi/2;pi/2)
 581   //We have 10 bins
 582   Int_t irp=Int_t(10.*fRPfull/TMath::Pi());
 583   if(irp>9)irp=9 ;
 584
 585   if(!fPHOSEvents[zvtx][fCenBin][irp])
 586     fPHOSEvents[zvtx][fCenBin][irp]=new TList() ;
 587   TList * prevPHOS = fPHOSEvents[zvtx][fCenBin][irp] ;
 588
 589   // Get PHOS rotation matrices from ESD and set them to the PHOS geometry
 590   if(fEventCounter == 0) {
 591     for(Int_t mod=0; mod<5; mod++) {
 592       const TGeoHMatrix* m =header->GetPHOSMatrix(mod) ;
 593       fPHOSGeo->SetMisalMatrix(m,mod) ;
 594       Printf("PHOS geo matrix for module # %d is set: %p\n", mod,m);
 595     }
 596     fEventCounter++ ;
 597   }
 598
 599   ProcessMC() ;
 600
 601   if(fPHOSEvent1){
 602     fPHOSEvent1->Clear() ;
 603     fPHOSEvent2->Clear() ;
 604   }
 605   else{
 606     fPHOSEvent1 = new TClonesArray("AliCaloPhoton",200) ;
 607     fPHOSEvent2 = new TClonesArray("AliCaloPhoton",200) ;
 608   }
 609
 610   TClonesArray * clustersEmb = (TClonesArray*)event->FindListObject("EmbeddedCaloClusters") ;
 611   AliAODCaloCells * cellsEmb = (AliAODCaloCells *)event->FindListObject("EmbeddedPHOScells") ;
 612   TClonesArray * clustersOld = event->GetCaloClusters() ;
 613   AliAODCaloCells * cellsOld = event->GetPHOSCells() ;
 614   TVector3 vertex(vtx0);
 615   char key[55] ;
 616   //Before Embedding
 617   Int_t multClustOld = clustersOld->GetEntriesFast();
 618   Int_t multClustEmb = clustersEmb->GetEntriesFast();
 619   Int_t inPHOSold=0 ;
 620   Int_t inPHOSemb=0 ;
 621   for (Int_t i=0; i<multClustOld; i++) {
 622     AliAODCaloCluster *clu = (AliAODCaloCluster*)clustersOld->At(i);
 623     if ( !clu->IsPHOS() || clu->E()<0.3) continue;
 624
 625     Bool_t survive=kFALSE ;
 626     for(Int_t ii=0;(ii<multClustEmb)&&(!survive);ii++){
 627        AliAODCaloCluster *clu2 = (AliAODCaloCluster*)clustersEmb->At(ii);
 628        survive=IsSameCluster(clu,clu2);
 629     }
 630
 631
 632     Float_t  position[3];
 633     clu->GetPosition(position);
 634     TVector3 global(position) ;
 635     Int_t relId[4] ;
 636     fPHOSGeo->GlobalPos2RelId(global,relId) ;
 637     Int_t mod  = relId[0] ;
 638     Int_t cellX = relId[2];
 639     Int_t cellZ = relId[3] ;
 640     if ( !IsGoodChannel("PHOS",mod,cellX,cellZ) )
 641       continue ;
 642     if(clu->GetNCells()<3)
 643       continue ;
 644
 645     snprintf(key,55,"hCluM%d",mod) ;
 646     FillHistogram(key,cellX,cellZ,1.);
 647
 648     TLorentzVector pv1 ;
 649     clu->GetMomentum(pv1 ,vtx0);
 650
 651     if(inPHOSold>=fPHOSEvent1->GetSize()){
 652       fPHOSEvent1->Expand(inPHOSold+50) ;
 653     }
 654     AliCaloPhoton * ph = new((*fPHOSEvent1)[inPHOSold]) AliCaloPhoton(pv1.X(),pv1.Py(),pv1.Z(),pv1.E()) ;
 655     ph->SetModule(mod) ;
 656     AliPHOSAodCluster cluPHOS1(*clu);
 657     cluPHOS1.Recalibrate(fPHOSCalibData,cellsOld); // modify the cell energies
 658     Double_t ecore=CoreEnergy(&cluPHOS1) ;
 659     pv1*= ecore/pv1.E() ;
 660     ph->SetMomV2(&pv1) ;
 661     ph->SetNCells(clu->GetNCells());
 662     ph->SetDispBit(TestLambda(clu->E(),clu->GetM20(),clu->GetM02())) ;
 663     ph->SetDisp2Bit(TestLambda2(clu->E(),clu->GetM20(),clu->GetM02())) ;
 664     ph->SetCPVBit(clu->GetEmcCpvDistance()>2.) ;
 665     ph->SetCPV2Bit(clu->GetEmcCpvDistance()>4.) ;
 666     if(!survive) //this cluster found in list after embedding, skipping it
 667       ph->SetTagged(1) ;
 668     ph->SetPhoton(clu->GetNExMax()<2); // Remember, if it is unfolded
 669
 670     if(!survive){
 671       Double_t distBC=clu->GetDistanceToBadChannel();
 672       if(distBC>2.)
 673         FillHistogram(Form("hPhotAll_DistBad2_cen%d",fCenBin),ph->Pt(),-1.) ;
 674         if(distBC>4.)
 675           FillHistogram(Form("hPhotAll_DistBad4_cen%d",fCenBin),ph->Pt(),-1.) ;
 676           if(distBC>6.)
 677             FillHistogram(Form("hPhotAll_DistBad6_cen%d",fCenBin),ph->Pt(),-1.) ;
 678     }
 679     inPHOSold++ ;
 680   }
 681
 682   for (Int_t i=0; i<multClustEmb; i++) {
 683     AliAODCaloCluster *clu = (AliAODCaloCluster*)clustersEmb->At(i);
 684     if ( !clu->IsPHOS() || clu->E()<0.3) continue;
 685
 686     Bool_t survive=kFALSE ;
 687     for(Int_t ii=0;(ii<multClustOld)&&(!survive);ii++){
 688        AliAODCaloCluster *clu2 = (AliAODCaloCluster*)clustersOld->At(ii);
 689        survive=IsSameCluster(clu,clu2);
 690     }
 691
 692     Float_t  position[3];
 693     clu->GetPosition(position);
 694     TVector3 global(position) ;
 695     Int_t relId[4] ;
 696     fPHOSGeo->GlobalPos2RelId(global,relId) ;
 697     Int_t mod  = relId[0] ;
 698     Int_t cellX = relId[2];
 699     Int_t cellZ = relId[3] ;
 700     if ( !IsGoodChannel("PHOS",mod,cellX,cellZ) )
 701       continue ;
 702     if(clu->GetNCells()<3)
 703       continue ;
 704
 705     snprintf(key,55,"hCluM%d",mod) ;
 706     FillHistogram(key,cellX,cellZ,1.);
 707
 708     TLorentzVector pv1 ;
 709     clu->GetMomentum(pv1 ,vtx0);
 710
 711     if(inPHOSemb>=fPHOSEvent2->GetSize()){
 712       fPHOSEvent2->Expand(inPHOSemb+50) ;
 713     }
 714     AliCaloPhoton * ph = new((*fPHOSEvent2)[inPHOSemb]) AliCaloPhoton(pv1.X(),pv1.Py(),pv1.Z(),pv1.E()) ;
 715     ph->SetModule(mod) ;
 716     AliPHOSAodCluster cluPHOS1(*clu);
 717     cluPHOS1.Recalibrate(fPHOSCalibData,cellsEmb); // modify the cell energies
 718     Double_t ecore=CoreEnergy(&cluPHOS1) ;
 719     pv1*= ecore/pv1.E() ;
 720     ph->SetMomV2(&pv1) ;
 721     ph->SetNCells(clu->GetNCells());
 722     ph->SetDispBit(TestLambda(clu->E(),clu->GetM20(),clu->GetM02())) ;
 723     ph->SetDisp2Bit(TestLambda2(clu->E(),clu->GetM20(),clu->GetM02())) ;
 724     ph->SetCPVBit(clu->GetEmcCpvDistance()>2.) ;
 725     ph->SetCPV2Bit(clu->GetEmcCpvDistance()>4.) ;
 726     if(!survive) //this cluster found in list after embedding, skipping it
 727       ph->SetTagged(1) ;
 728     ph->SetPhoton(clu->GetNExMax()<2); // Remember, if it is unfolded
 729
 730     if(!survive){
 731       Double_t distBC=clu->GetDistanceToBadChannel();
 732       if(distBC>2.)
 733         FillHistogram(Form("hPhotAll_DistBad2_cen%d",fCenBin),ph->Pt()) ;
 734         if(distBC>4.)
 735           FillHistogram(Form("hPhotAll_DistBad4_cen%d",fCenBin),ph->Pt()) ;
 736           if(distBC>6.)
 737             FillHistogram(Form("hPhotAll_DistBad6_cen%d",fCenBin),ph->Pt()) ;
 738     }
 739
 740     inPHOSemb++ ;
 741   }
 742
 743   //Single photon
 744   for (Int_t i1=0; i1<inPHOSold; i1++) {
 745     AliCaloPhoton * ph1=(AliCaloPhoton*)fPHOSEvent1->At(i1) ;
 746     if(!ph1->IsTagged())
 747       continue ;
 748     snprintf(key,55,"hPhotAll_cen%d",fCenBin) ;
 749     FillHistogram(key,ph1->Pt(),-1.) ;
 750     snprintf(key,55,"hNegPhotAll_cen%d",fCenBin) ;
 751     FillHistogram(key,ph1->Pt(),-1.) ;
 752     snprintf(key,55,"hPhotAllcore_cen%d",fCenBin) ;
 753     FillHistogram(key,ph1->GetMomV2()->Pt(),-1.) ;
 754     snprintf(key,55,"hNegPhotAllcore_cen%d",fCenBin) ;
 755     FillHistogram(key,ph1->GetMomV2()->Pt(),-1.) ;
 756     if(ph1->IsPhoton()){
 757       snprintf(key,55,"hPhotAllwou_cen%d",fCenBin) ;
 758       FillHistogram(key,ph1->Pt(),-1.) ;
 759       snprintf(key,55,"hNegPhotAllwou_cen%d",fCenBin) ;
 760       FillHistogram(key,ph1->Pt(),-1.) ;
 761     }
 762     if(ph1->IsCPVOK() ){
 763       snprintf(key,55,"hPhotCPV_cen%d",fCenBin) ;
 764       FillHistogram(key,ph1->Pt(),-1.) ;
 765       snprintf(key,55,"hNegPhotCPV_cen%d",fCenBin) ;
 766       FillHistogram(key,ph1->Pt(),-1.) ;
 767       snprintf(key,55,"hPhotCPVcore_cen%d",fCenBin) ;
 768       FillHistogram(key,ph1->GetMomV2()->Pt(),-1.) ;
 769       snprintf(key,55,"hNegPhotCPVcore_cen%d",fCenBin) ;
 770       FillHistogram(key,ph1->GetMomV2()->Pt(),-1.) ;
 771     }
 772     if(ph1->IsCPV2OK() ){
 773       snprintf(key,55,"hPhotCPV2_cen%d",fCenBin) ;
 774       FillHistogram(key,ph1->Pt(),-1.) ;
 775       snprintf(key,55,"hNegPhotCPV2_cen%d",fCenBin) ;
 776       FillHistogram(key,ph1->Pt(),-1.) ;
 777     }
 778     if(ph1->IsDisp2OK()){
 779       snprintf(key,55,"hPhotDisp2_cen%d",fCenBin) ;
 780       FillHistogram(key,ph1->Pt(),-1.) ;
 781       snprintf(key,55,"hNegPhotDisp2_cen%d",fCenBin) ;
 782       FillHistogram(key,ph1->Pt(),-1.) ;
 783     }
 784     if(ph1->IsDispOK()){
 785       snprintf(key,55,"hPhotDisp_cen%d",fCenBin) ;
 786       FillHistogram(key,ph1->Pt(),-1.) ;
 787       snprintf(key,55,"hNegPhotDisp_cen%d",fCenBin) ;
 788       FillHistogram(key,ph1->Pt(),-1.) ;
 789       snprintf(key,55,"hPhotDispcore_cen%d",fCenBin) ;
 790       FillHistogram(key,ph1->GetMomV2()->Pt(),-1.) ;
 791       snprintf(key,55,"hNegPhotDispcore_cen%d",fCenBin) ;
 792       FillHistogram(key,ph1->GetMomV2()->Pt(),-1.) ;
 793       if(ph1->IsCPVOK()){
 794         snprintf(key,55,"hPhotBoth_cen%d",fCenBin) ;
 795         FillHistogram(key,ph1->Pt(),-1.) ;
 796         snprintf(key,55,"hNegPhotBoth_cen%d",fCenBin) ;
 797         FillHistogram(key,ph1->Pt(),-1.) ;
 798         snprintf(key,55,"hPhotBothcore_cen%d",fCenBin) ;
 799         FillHistogram(key,ph1->GetMomV2()->Pt(),-1.) ;
 800         snprintf(key,55,"hNegPhotBothcore_cen%d",fCenBin) ;
 801         FillHistogram(key,ph1->GetMomV2()->Pt(),-1.) ;
 802       }
 803     } // end of loop i2
 804   } // end of loop i1
 805
 806   for (Int_t i1=0; i1<inPHOSemb; i1++) {
 807     AliCaloPhoton * ph1=(AliCaloPhoton*)fPHOSEvent2->At(i1) ;
 808     if(!ph1->IsTagged())
 809       continue ;
 810     snprintf(key,55,"hPhotAll_cen%d",fCenBin) ;
 811     FillHistogram(key,ph1->Pt(),1.) ;
 812     snprintf(key,55,"hPhotAllcore_cen%d",fCenBin) ;
 813     FillHistogram(key,ph1->GetMomV2()->Pt(),1.) ;
 814     if(ph1->IsPhoton()){
 815       snprintf(key,55,"hPhotAllwou_cen%d",fCenBin) ;
 816       FillHistogram(key,ph1->Pt(),1.) ;
 817     }
 818     if(ph1->IsCPVOK() ){
 819       snprintf(key,55,"hPhotCPV_cen%d",fCenBin) ;
 820       FillHistogram(key,ph1->Pt(),1.) ;
 821       snprintf(key,55,"hPhotCPVcore_cen%d",fCenBin) ;
 822       FillHistogram(key,ph1->GetMomV2()->Pt(),1.) ;
 823     }
 824     if(ph1->IsCPV2OK() ){
 825       snprintf(key,55,"hPhotCPV2_cen%d",fCenBin) ;
 826       FillHistogram(key,ph1->Pt(),1.) ;
 827     }
 828     if(ph1->IsDisp2OK()){
 829       snprintf(key,55,"hPhotDisp2_cen%d",fCenBin) ;
 830       FillHistogram(key,ph1->Pt(),1.) ;
 831     }
 832     if(ph1->IsDispOK()){
 833       snprintf(key,55,"hPhotDisp_cen%d",fCenBin) ;
 834       FillHistogram(key,ph1->Pt(),1.) ;
 835       snprintf(key,55,"hPhotDispcore_cen%d",fCenBin) ;
 836       FillHistogram(key,ph1->GetMomV2()->Pt(),1.) ;
 837       if(ph1->IsCPVOK()){
 838         snprintf(key,55,"hPhotBoth_cen%d",fCenBin) ;
 839         FillHistogram(key,ph1->Pt(),1.) ;
 840         snprintf(key,55,"hPhotBothcore_cen%d",fCenBin) ;
 841         FillHistogram(key,ph1->GetMomV2()->Pt(),1.) ;
 842       }
 843     } // end of loop i2
 844   } // end of loop i1
 845
 846
 847
 848   // Fill Real disribution:
 849   // Disappeared clusters enter with negative contribution
 850   // In addition fill control histogam with Real before embedding
 851   for (Int_t i1=0; i1<inPHOSold-1; i1++) {
 852     AliCaloPhoton * ph1=(AliCaloPhoton*)fPHOSEvent1->At(i1) ;
 853     for (Int_t i2=i1+1; i2<inPHOSold; i2++) {
 854       AliCaloPhoton * ph2=(AliCaloPhoton*)fPHOSEvent1->At(i2) ;
 855
 856       TLorentzVector p12  = *ph1  + *ph2;
 857       TLorentzVector pv12 = *(ph1->GetMomV2()) + *(ph2->GetMomV2());
 858       Double_t a=TMath::Abs((ph1->E()-ph2->E())/(ph1->E()+ph2->E())) ;
 859
 860       //Fill Controll histogram: Real before embedding
 861       snprintf(key,55,"hOldMassPtAll_cen%d",fCenBin) ;
 862       FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 863       snprintf(key,55,"hOldMassPtAllcore_cen%d",fCenBin) ;
 864       FillHistogram(key,pv12.M() ,pv12.Pt(),-1.) ;
 865       if(ph1->IsPhoton() && ph2->IsPhoton()){
 866         snprintf(key,55,"hOldMassPtAllwou_cen%d",fCenBin) ;
 867         FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 868       }
 869       if(ph1->IsCPVOK() && ph2->IsCPVOK()){
 870         snprintf(key,55,"hOldMassPtCPV_cen%d",fCenBin) ;
 871         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 872         snprintf(key,55,"hOldMassPtCPV_cen%d",fCenBin) ;
 873         FillHistogram(key,pv12.M(), pv12.Pt(),-1) ;
 874       }
 875       if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
 876         snprintf(key,55,"hOldMassPtCPV2_cen%d",fCenBin) ;
 877         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 878       }
 879       if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
 880         snprintf(key,55,"hOldMassPtDisp2_cen%d",fCenBin) ;
 881         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 882       }
 883       if(ph1->IsDispOK() && ph2->IsDispOK()){
 884         snprintf(key,55,"hOldMassPtDisp_cen%d",fCenBin) ;
 885         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 886         snprintf(key,55,"hOldMassPtDispcore_cen%d",fCenBin) ;
 887         FillHistogram(key,pv12.M() ,pv12.Pt(),-1) ;
 888         if(ph1->IsCPVOK() && ph2->IsCPVOK()){
 889           snprintf(key,55,"hOldMassPtBoth_cen%d",fCenBin) ;
 890           FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 891           snprintf(key,55,"hOldMassPtBothcore_cen%d",fCenBin) ;
 892           FillHistogram(key,pv12.M() ,pv12.Pt(),-1) ;
 893         }
 894       }
 895
 896       //Now fill main histograms with negative contributions
 897       if(!(ph1->IsTagged() || ph2->IsTagged()) )
 898         continue ;
 899       snprintf(key,55,"hMassPtAll_cen%d",fCenBin) ;
 900       FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 901       snprintf(key,55,"hMassPtAllcore_cen%d",fCenBin) ;
 902       FillHistogram(key,pv12.M() ,pv12.Pt(),-1.) ;
 903       if(ph1->IsPhoton() && ph2->IsPhoton()){
 904         snprintf(key,55,"hMassPtAllwou_cen%d",fCenBin) ;
 905         FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 906       }
 907       if(a<0.9){
 908         snprintf(key,55,"hMassPtAll_a09_cen%d",fCenBin) ;
 909         FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 910         if(a<0.8){
 911           snprintf(key,55,"hMassPtAll_a08_cen%d",fCenBin) ;
 912           FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 913           if(a<0.7){
 914             snprintf(key,55,"hMassPtAll_a07_cen%d",fCenBin) ;
 915             FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 916           }
 917         }
 918       }
 919       snprintf(key,55,"hNegMassPtAll_cen%d",fCenBin) ;
 920       FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 921       snprintf(key,55,"hNegMassPtAllcore_cen%d",fCenBin) ;
 922       FillHistogram(key,pv12.M() ,pv12.Pt(),-1.) ;
 923       if(ph1->IsPhoton() && ph2->IsPhoton()){
 924         snprintf(key,55,"hNegMassPtAllwou_cen%d",fCenBin) ;
 925         FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 926       }
 927
 928       if(ph1->IsCPVOK() && ph2->IsCPVOK()){
 929         snprintf(key,55,"hMassPtCPV_cen%d",fCenBin) ;
 930         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 931         snprintf(key,55,"hMassPtCPVcore_cen%d",fCenBin) ;
 932         FillHistogram(key,pv12.M() ,pv12.Pt(),-1) ;
 933         if(a<0.9){
 934           snprintf(key,55,"hMassPtCPV_a09_cen%d",fCenBin) ;
 935           FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 936           if(a<0.8){
 937             snprintf(key,55,"hMassPtCPV_a08_cen%d",fCenBin) ;
 938             FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 939             if(a<0.7){
 940               snprintf(key,55,"hMassPtCPV_a07_cen%d",fCenBin) ;
 941               FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 942             }
 943           }
 944         }
 945         snprintf(key,55,"hNegMassPtCPV_cen%d",fCenBin) ;
 946         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 947         snprintf(key,55,"hNegMassPtCPVcore_cen%d",fCenBin) ;
 948         FillHistogram(key,pv12.M() ,pv12.Pt(),-1) ;
 949       }
 950       if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
 951         snprintf(key,55,"hMassPtCPV2_cen%d",fCenBin) ;
 952         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 953         if(a<0.9){
 954           snprintf(key,55,"hMassPtCPV2_a09_cen%d",fCenBin) ;
 955           FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 956           if(a<0.8){
 957             snprintf(key,55,"hMassPtCPV2_a08_cen%d",fCenBin) ;
 958             FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 959             if(a<0.7){
 960               snprintf(key,55,"hMassPtCPV2_a07_cen%d",fCenBin) ;
 961               FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 962             }
 963           }
 964         }
 965         snprintf(key,55,"hNegMassPtCPV2_cen%d",fCenBin) ;
 966         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 967       }
 968
 969       if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
 970         snprintf(key,55,"hMassPtDisp2_cen%d",fCenBin) ;
 971         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 972       }
 973       if(ph1->IsDispOK() && ph2->IsDispOK()){
 974         snprintf(key,55,"hMassPtDisp_cen%d",fCenBin) ;
 975         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 976         snprintf(key,55,"hMassPtDispcore_cen%d",fCenBin) ;
 977         FillHistogram(key,pv12.M() ,pv12.Pt(),-1) ;
 978         if(a<0.9){
 979           snprintf(key,55,"hMassPtDisp_a09_cen%d",fCenBin) ;
 980           FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 981           if(a<0.8){
 982             snprintf(key,55,"hMassPtDisp_a08_cen%d",fCenBin) ;
 983             FillHistogram(key,p12.M() ,p12.Pt()) ;
 984             if(a<0.7){
 985               snprintf(key,55,"hMassPtDisp_a07_cen%d",fCenBin) ;
 986               FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
 987             }
 988           }
 989         }
 990         snprintf(key,55,"hNegMassPtDisp_cen%d",fCenBin) ;
 991         FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 992         snprintf(key,55,"hNegMassPtDispcore_cen%d",fCenBin) ;
 993         FillHistogram(key,pv12.M() ,pv12.Pt(),-1) ;
 994         if(ph1->IsCPVOK() && ph2->IsCPVOK()){
 995           snprintf(key,55,"hMassPtBoth_cen%d",fCenBin) ;
 996           FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
 997           snprintf(key,55,"hMassPtBothcore_cen%d",fCenBin) ;
 998           FillHistogram(key,pv12.M() ,pv12.Pt(),-1) ;
 999           if(a<0.9){
1000             snprintf(key,55,"hMassPtBoth_a09_cen%d",fCenBin) ;
1001             FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
1002             if(a<0.8){
1003               snprintf(key,55,"hMassPtBoth_a08_cen%d",fCenBin) ;
1004               FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
1005               if(a<0.7){
1006                 snprintf(key,55,"hMassPtBoth_a07_cen%d",fCenBin) ;
1007                 FillHistogram(key,p12.M() ,p12.Pt(),-1.) ;
1008              }
1009             }
1010           }
1011           snprintf(key,55,"hNegMassPtBoth_cen%d",fCenBin) ;
1012           FillHistogram(key,p12.M() ,p12.Pt(),-1) ;
1013           snprintf(key,55,"hNegMassPtBothcore_cen%d",fCenBin) ;
1014           FillHistogram(key,pv12.M() ,pv12.Pt(),-1) ;
1015         }
1016       }
1017     } // end of loop i2
1018   } // end of loop i1
1019
1020
1021   // Further fill Real disribution
1022   // now with positive contribution from new clusters
1023   // ass well fill controll histogram
1024   for (Int_t i1=0; i1<inPHOSemb-1; i1++) {
1025     AliCaloPhoton * ph1=(AliCaloPhoton*)fPHOSEvent2->At(i1) ;
1026     for (Int_t i2=i1+1; i2<inPHOSemb; i2++) {
1027       AliCaloPhoton * ph2=(AliCaloPhoton*)fPHOSEvent2->At(i2) ;
1028
1029       TLorentzVector p12  = *ph1  + *ph2;
1030       TLorentzVector pv12 = *(ph1->GetMomV2()) + *(ph2->GetMomV2());
1031       Double_t a=TMath::Abs((ph1->E()-ph2->E())/(ph1->E()+ph2->E())) ;
1032
1033       // Controll histogram: Real after embedding
1034       snprintf(key,55,"hNewMassPtAll_cen%d",fCenBin) ;
1035       FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1036       snprintf(key,55,"hNewMassPtAllcore_cen%d",fCenBin) ;
1037       FillHistogram(key,pv12.M() ,pv12.Pt(),1.) ;
1038       if(ph1->IsPhoton() && ph2->IsPhoton()){
1039         snprintf(key,55,"hNewMassPtAllwou_cen%d",fCenBin) ;
1040         FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1041       }
1042       if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1043         snprintf(key,55,"hNewMassPtCPV_cen%d",fCenBin) ;
1044         FillHistogram(key,p12.M() ,p12.Pt(),1) ;
1045         snprintf(key,55,"hNewMassPtCPVcore_cen%d",fCenBin) ;
1046         FillHistogram(key,pv12.M() ,pv12.Pt(),1) ;
1047       }
1048       if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
1049         snprintf(key,55,"hNewMassPtCPV2_cen%d",fCenBin) ;
1050         FillHistogram(key,p12.M() ,p12.Pt(),1) ;
1051       }
1052       if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
1053         snprintf(key,55,"hNewMassPtDisp2_cen%d",fCenBin) ;
1054         FillHistogram(key,p12.M() ,p12.Pt(),1) ;
1055       }
1056       if(ph1->IsDispOK() && ph2->IsDispOK()){
1057         snprintf(key,55,"hNewMassPtDisp_cen%d",fCenBin) ;
1058         FillHistogram(key,p12.M() ,p12.Pt(),1) ;
1059         snprintf(key,55,"hNewMassPtDispcore_cen%d",fCenBin) ;
1060         FillHistogram(key,pv12.M() ,pv12.Pt(),1) ;
1061         if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1062           snprintf(key,55,"hNewMassPtBoth_cen%d",fCenBin) ;
1063           FillHistogram(key,p12.M() ,p12.Pt(),1) ;
1064           snprintf(key,55,"hNewMassPtBothcore_cen%d",fCenBin) ;
1065           FillHistogram(key,pv12.M() ,pv12.Pt(),1) ;
1066         }
1067       }
1068
1069       //Now fill main histogamm
1070       //new clusters with positive contribution
1071       if(!(ph1->IsTagged() || ph2->IsTagged()) )
1072         continue ;
1073       snprintf(key,55,"hMassPtAll_cen%d",fCenBin) ;
1074       FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1075       snprintf(key,55,"hMassPtAllcore_cen%d",fCenBin) ;
1076       FillHistogram(key,pv12.M() ,pv12.Pt(),1.) ;
1077       if(ph1->IsPhoton() && ph2->IsPhoton()){
1078         snprintf(key,55,"hMassPtAllwou_cen%d",fCenBin) ;
1079         FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1080       }
1081       if(a<0.9){
1082         snprintf(key,55,"hMassPtAll_a09_cen%d",fCenBin) ;
1083         FillHistogram(key,p12.M() ,p12.Pt()) ;
1084         if(a<0.8){
1085           snprintf(key,55,"hMassPtAll_a08_cen%d",fCenBin) ;
1086           FillHistogram(key,p12.M() ,p12.Pt()) ;
1087           if(a<0.7){
1088             snprintf(key,55,"hMassPtAll_a07_cen%d",fCenBin) ;
1089             FillHistogram(key,p12.M() ,p12.Pt()) ;
1090           }
1091         }
1092       }
1093       if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1094         snprintf(key,55,"hMassPtCPV_cen%d",fCenBin) ;
1095         FillHistogram(key,p12.M() ,p12.Pt(),1) ;
1096         snprintf(key,55,"hMassPtCPVcore_cen%d",fCenBin) ;
1097         FillHistogram(key,pv12.M() ,pv12.Pt(),1) ;
1098         if(a<0.9){
1099           snprintf(key,55,"hMassPtCPV_a09_cen%d",fCenBin) ;
1100           FillHistogram(key,p12.M() ,p12.Pt()) ;
1101           if(a<0.8){
1102             snprintf(key,55,"hMassPtCPV_a08_cen%d",fCenBin) ;
1103             FillHistogram(key,p12.M() ,p12.Pt()) ;
1104             if(a<0.7){
1105               snprintf(key,55,"hMassPtCPV_a07_cen%d",fCenBin) ;
1106               FillHistogram(key,p12.M() ,p12.Pt()) ;
1107             }
1108           }
1109         }
1110       }
1111       if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
1112         snprintf(key,55,"hMassPtCPV2_cen%d",fCenBin) ;
1113         FillHistogram(key,p12.M() ,p12.Pt(),1) ;
1114          if(a<0.9){
1115           snprintf(key,55,"hMassPtCPV2_a09_cen%d",fCenBin) ;
1116           FillHistogram(key,p12.M() ,p12.Pt()) ;
1117           if(a<0.8){
1118             snprintf(key,55,"hMassPtCPV2_a08_cen%d",fCenBin) ;
1119             FillHistogram(key,p12.M() ,p12.Pt()) ;
1120             if(a<0.7){
1121               snprintf(key,55,"hMassPtCPV2_a07_cen%d",fCenBin) ;
1122               FillHistogram(key,p12.M() ,p12.Pt()) ;
1123             }
1124           }
1125         }
1126       }
1127       if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
1128         snprintf(key,55,"hMassPtDisp2_cen%d",fCenBin) ;
1129         FillHistogram(key,p12.M() ,p12.Pt(),1) ;
1130       }
1131       if(ph1->IsDispOK() && ph2->IsDispOK()){
1132         snprintf(key,55,"hMassPtDisp_cen%d",fCenBin) ;
1133         FillHistogram(key,p12.M() ,p12.Pt(),1) ;
1134         snprintf(key,55,"hMassPtDispcore_cen%d",fCenBin) ;
1135         FillHistogram(key,pv12.M() ,pv12.Pt(),1) ;
1136         if(a<0.9){
1137           snprintf(key,55,"hMassPtDisp_a09_cen%d",fCenBin) ;
1138           FillHistogram(key,p12.M() ,p12.Pt()) ;
1139           if(a<0.8){
1140             snprintf(key,55,"hMassPtDisp_a08_cen%d",fCenBin) ;
1141             FillHistogram(key,p12.M() ,p12.Pt()) ;
1142             if(a<0.7){
1143               snprintf(key,55,"hMassPtDisp_a07_cen%d",fCenBin) ;
1144               FillHistogram(key,p12.M() ,p12.Pt()) ;
1145             }
1146           }
1147         }
1148
1149         if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1150           snprintf(key,55,"hMassPtBoth_cen%d",fCenBin) ;
1151           FillHistogram(key,p12.M() ,p12.Pt(),1) ;
1152           snprintf(key,55,"hMassPtBothcore_cen%d",fCenBin) ;
1153           FillHistogram(key,pv12.M() ,pv12.Pt(),1) ;
1154           if(a<0.9){
1155             snprintf(key,55,"hMassPtBoth_a09_cen%d",fCenBin) ;
1156             FillHistogram(key,p12.M() ,p12.Pt()) ;
1157             if(a<0.8){
1158               snprintf(key,55,"hMassPtBoth_a08_cen%d",fCenBin) ;
1159               FillHistogram(key,p12.M() ,p12.Pt()) ;
1160               if(a<0.7){
1161                 snprintf(key,55,"hMassPtBoth_a07_cen%d",fCenBin) ;
1162                 FillHistogram(key,p12.M() ,p12.Pt()) ;
1163               }
1164             }
1165           }
1166         }
1167       }
1168     } // end of loop i2
1169   } // end of loop i1
1170
1171
1172   //now mixed, does not really matter old or new list of clusters
1173   for (Int_t i1=0; i1<inPHOSemb; i1++) {
1174     AliCaloPhoton * ph1=(AliCaloPhoton*)fPHOSEvent2->At(i1) ;
1175     for(Int_t ev=0; ev<prevPHOS->GetSize();ev++){
1176       TClonesArray * mixPHOS = static_cast<TClonesArray*>(prevPHOS->At(ev)) ;
1177       for(Int_t i2=0; i2<mixPHOS->GetEntriesFast();i2++){
1178         AliCaloPhoton * ph2=(AliCaloPhoton*)mixPHOS->At(i2) ;
1179         TLorentzVector p12  = *ph1  + *ph2;
1180         TLorentzVector pv12 = *(ph1->GetMomV2()) + *(ph2->GetMomV2());
1181         Double_t a=TMath::Abs((ph1->E()-ph2->E())/(ph1->E()+ph2->E())) ;
1182
1183         snprintf(key,55,"hMiMassPtAll_cen%d",fCenBin) ;
1184         FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1185         snprintf(key,55,"hMiMassPtAllcore_cen%d",fCenBin) ;
1186         FillHistogram(key,pv12.M() ,pv12.Pt(),1.) ;
1187         if(ph1->IsPhoton() && ph2->IsPhoton()){
1188           snprintf(key,55,"hMiMassPtAllwou_cen%d",fCenBin) ;
1189           FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1190         }
1191         if(a<0.9){
1192           snprintf(key,55,"hMiMassPtAll_a09_cen%d",fCenBin) ;
1193           FillHistogram(key,p12.M() ,p12.Pt()) ;
1194           if(a<0.8){
1195             snprintf(key,55,"hMiMassPtAll_a08_cen%d",fCenBin) ;
1196             FillHistogram(key,p12.M() ,p12.Pt()) ;
1197             if(a<0.7){
1198               snprintf(key,55,"hMiMassPtAll_a07_cen%d",fCenBin) ;
1199               FillHistogram(key,p12.M() ,p12.Pt()) ;
1200             }
1201           }
1202         }
1203         if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1204           snprintf(key,55,"hMiMassPtCPV_cen%d",fCenBin) ;
1205           FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1206           snprintf(key,55,"hMiMassPtCPVcore_cen%d",fCenBin) ;
1207           FillHistogram(key,pv12.M() ,pv12.Pt(),1.) ;
1208           if(a<0.9){
1209             snprintf(key,55,"hMiMassPtCPV_a09_cen%d",fCenBin) ;
1210             FillHistogram(key,p12.M() ,p12.Pt()) ;
1211             if(a<0.8){
1212               snprintf(key,55,"hMiMassPtCPV_a08_cen%d",fCenBin) ;
1213               FillHistogram(key,p12.M() ,p12.Pt()) ;
1214               if(a<0.7){
1215                 snprintf(key,55,"hMiMassPtCPV_a07_cen%d",fCenBin) ;
1216                 FillHistogram(key,p12.M() ,p12.Pt()) ;
1217               }
1218             }
1219           }
1220         }
1221         if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
1222           snprintf(key,55,"hMiMassPtCPV2_cen%d",fCenBin) ;
1223           FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1224           if(a<0.9){
1225             snprintf(key,55,"hMiMassPtCPV2_a09_cen%d",fCenBin) ;
1226             FillHistogram(key,p12.M() ,p12.Pt()) ;
1227             if(a<0.8){
1228               snprintf(key,55,"hMiMassPtCPV2_a08_cen%d",fCenBin) ;
1229               FillHistogram(key,p12.M() ,p12.Pt()) ;
1230               if(a<0.7){
1231                 snprintf(key,55,"hMiMassPtCPV2_a07_cen%d",fCenBin) ;
1232                 FillHistogram(key,p12.M() ,p12.Pt()) ;
1233               }
1234             }
1235           }
1236         }
1237         if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
1238           snprintf(key,55,"hMiMassPtDisp2_cen%d",fCenBin) ;
1239           FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1240         }
1241         if(ph1->IsDispOK() && ph2->IsDispOK()){
1242           snprintf(key,55,"hMiMassPtDisp_cen%d",fCenBin) ;
1243           FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1244           snprintf(key,55,"hMiMassPtDispcore_cen%d",fCenBin) ;
1245           FillHistogram(key,pv12.M() ,pv12.Pt(),1.) ;
1246           if(a<0.9){
1247             snprintf(key,55,"hMiMassPtDisp_a09_cen%d",fCenBin) ;
1248             FillHistogram(key,p12.M() ,p12.Pt()) ;
1249             if(a<0.8){
1250               snprintf(key,55,"hMiMassPtDisp_a08_cen%d",fCenBin) ;
1251               FillHistogram(key,p12.M() ,p12.Pt()) ;
1252               if(a<0.7){
1253                 snprintf(key,55,"hMiMassPtDisp_a07_cen%d",fCenBin) ;
1254                 FillHistogram(key,p12.M() ,p12.Pt()) ;
1255               }
1256             }
1257           }
1258
1259           if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1260             snprintf(key,55,"hMiMassPtBoth_cen%d",fCenBin) ;
1261             FillHistogram(key,p12.M() ,p12.Pt(),1.) ;
1262             snprintf(key,55,"hMiMassPtBothcore_cen%d",fCenBin) ;
1263             FillHistogram(key,pv12.M() ,pv12.Pt(),1.) ;
1264             if(a<0.9){
1265               snprintf(key,55,"hMiMassPtBoth_a09_cen%d",fCenBin) ;
1266               FillHistogram(key,p12.M() ,p12.Pt()) ;
1267               if(a<0.8){
1268                 snprintf(key,55,"hMiMassPtBoth_a08_cen%d",fCenBin) ;
1269                 FillHistogram(key,p12.M() ,p12.Pt()) ;
1270                 if(a<0.7){
1271                   snprintf(key,55,"hMiMassPtBoth_a07_cen%d",fCenBin) ;
1272                   FillHistogram(key,p12.M() ,p12.Pt()) ;
1273                }
1274               }
1275             }
1276           }
1277         }
1278       } // end of loop i2
1279     }
1280   } // end of loop i1
1281
1282
1283   //Now we either add current events to stack or remove
1284   //If no photons in current event - no need to add it to mixed
1285   if(fPHOSEvent2->GetEntriesFast()>0){
1286     prevPHOS->AddFirst(fPHOSEvent2) ;
1287     fPHOSEvent2=0;
1288     delete fPHOSEvent1;
1289     fPHOSEvent1=0;
1290     if(prevPHOS->GetSize()>100){//Remove redundant events
1291       TClonesArray * tmp = static_cast<TClonesArray*>(prevPHOS->Last()) ;
1292       prevPHOS->RemoveLast() ;
1293       delete tmp ;
1294     }
1295   }
1296   // Post output data.
1297   PostData(1, fOutputContainer);
1298   fEventCounter++;
1299 }
1300 //___________________________________________________________________________
1301 Bool_t AliAnalysisTaskPi0DiffEfficiency::IsSameCluster(AliAODCaloCluster * c1, AliAODCaloCluster * c2)const{
1302  //Compare clusters before and after embedding
1303  //clusters are the same if
1304  // - Energy changed less than 0.1%  (numerical accuracy in reconstruction)
1305  // - lists of digits are the same
1306
1307  if(c1->GetNCells() != c2->GetNCells())
1308    return kFALSE ;
1309
1310  if(TMath::Abs(c1->E()-c2->E())>0.01*c1->E())
1311    return kFALSE ;
1312
1313  UShort_t *list1 = c1->GetCellsAbsId() ;
1314  UShort_t *list2 = c2->GetCellsAbsId() ;
1315  for(Int_t i=0; i<c1->GetNCells(); i++){
1316   if(list1[i] != list2[i])
1317     return kFALSE ;
1318  }
1319  return kTRUE ;
1320
1321 }
1322
1323
1324