PWG4/PartCorrDep/AliAnaOmegaToPi0Gamma.cxx

   1 /**************************************************************************\r
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
   3  *                                                                        *\r
   4  * Author: The ALICE Off-line Project.                                    *\r
   5  * Contributors are mentioned in the code where appropriate.              *\r
   6  *                                                                        *\r
   7  * Permission to use, copy, modify and distribute this software and its   *\r
   8  * documentation strictly for non-commercial purposes is hereby granted   *\r
   9  * without fee, provided that the above copyright notice appears in all   *\r
  10  * copies and that both the copyright notice and this permission notice   *\r
  11  * appear in the supporting documentation. The authors make no claims     *\r
  12  * about the suitability of this software for any purpose. It is          *\r
  13  * provided "as is" without express or implied warranty.                  *\r
  14  **************************************************************************/\r
  15 /* $Id: $ */\r
  16 //_________________________________________________________________________\r
  17 // class to extract omega(782)->pi0+gamma->3gamma\r
  18 //\r
  19 //-- Author: Renzhuo Wan (IOPP-Wuhan, China)\r
  20 //_________________________________________________________________________\r
  21 \r
  22 // --- ROOT system\r
  23 class TROOT;\r
  24 \r
  25 // --- AliRoot system\r
  26 //class AliVEvent;\r
  27 // --- ROOT system ---\r
  28 #include "TH2F.h"\r
  29 #include "TLorentzVector.h"\r
  30 #include "TParticle.h"\r
  31 #include "TCanvas.h"\r
  32 #include "TFile.h"\r
  33 //---- AliRoot system ----\r
  34 #include "AliAnaOmegaToPi0Gamma.h"\r
  35 #include "AliCaloTrackReader.h"\r
  36 #include "AliCaloPID.h"\r
  37 #include "AliStack.h"\r
  38 #include "AliVEvent.h"\r
  39 #include "AliAODEvent.h"\r
  40 #include "AliAODMCParticle.h"\r
  41 ClassImp(AliAnaOmegaToPi0Gamma)\r
  42 \r
  43 //______________________________________________________________________________\r
  44 AliAnaOmegaToPi0Gamma::AliAnaOmegaToPi0Gamma() : AliAnaPartCorrBaseClass(),\r
  45 fInputAODGamma(0),fInputAODPi0(0), fInputAODGammaName(""),\r
  46 fEventsList(0), \r
  47 fVtxZCut(0), fCent(0), fRp(0), fBadChDist(0),\r
  48 fNVtxZBin(0), fNCentBin(0), fNRpBin(0), fNBadChDistBin(0), fNpid(0),\r
  49 fNmaxMixEv(0), fPi0Mass(0),\r
  50 fPi0MassWindow(0),\r
  51 fNbinsPt(0),fPtBegin(0),fPtEnd(0),\r
  52 fNbinsM(0),fMinM(0),fMaxM(0),\r
  53 fhEtalon(0),\r
  54 fRealOmega(0), fMixAOmega(0),\r
  55 fMixBOmega(0), fMixCOmega(0),\r
  56 fRealOmega1(0), fMixAOmega1(0),\r
  57 fMixBOmega1(0), fMixCOmega1(0),\r
  58 fRealOmega2(0), fMixAOmega2(0),\r
  59 fMixBOmega2(0), fMixCOmega2(0),\r
  60 fhOmegaPriPt(0)\r
  61 {\r
  62  //Default Ctor\r
  63  InitParameters();\r
  64 }\r
  65 \r
  66 //______________________________________________________________________________\r
  67 AliAnaOmegaToPi0Gamma::AliAnaOmegaToPi0Gamma(const AliAnaOmegaToPi0Gamma & ex) : AliAnaPartCorrBaseClass(ex),\r
  68 fInputAODGamma(new TClonesArray (*ex.fInputAODGamma)), \r
  69 fInputAODPi0(new TClonesArray (*ex.fInputAODPi0)),\r
  70 fInputAODGammaName(ex.fInputAODGammaName),\r
  71 fEventsList(ex.fEventsList), \r
  72 fVtxZCut(ex.fVtxZCut), fCent(ex.fCent), fRp(ex.fRp), fBadChDist(ex.fBadChDist),\r
  73 fNVtxZBin(ex.fNVtxZBin), fNCentBin(ex.fNCentBin), fNRpBin(ex.fNRpBin),\r
  74 fNBadChDistBin(ex.fNBadChDistBin),fNpid(ex.fNpid),\r
  75 fNmaxMixEv(ex.fNmaxMixEv), fPi0Mass(ex.fPi0Mass),\r
  76 fPi0MassWindow(ex.fPi0MassWindow),\r
  77 fNbinsPt(ex.fNbinsPt),fPtBegin(ex.fPtBegin),fPtEnd(ex.fPtEnd),\r
  78 fNbinsM(ex.fNbinsM),fMinM(ex.fMinM),fMaxM(ex.fMaxM),\r
  79 fhEtalon(ex.fhEtalon),\r
  80 fRealOmega(ex.fRealOmega), fMixAOmega(ex.fMixAOmega),\r
  81 fMixBOmega(ex.fMixBOmega), fMixCOmega(ex.fMixCOmega),\r
  82 fRealOmega1(ex.fRealOmega1), fMixAOmega1(ex.fMixAOmega1),\r
  83 fMixBOmega1(ex.fMixBOmega1), fMixCOmega1(ex.fMixCOmega1),\r
  84 fRealOmega2(ex.fRealOmega2), fMixAOmega2(ex.fMixAOmega2),\r
  85 fMixBOmega2(ex.fMixBOmega2), fMixCOmega2(ex.fMixCOmega2),\r
  86 fhOmegaPriPt(ex.fhOmegaPriPt)\r
  87 {\r
  88  // cpy ctor\r
  89  //Do not need it\r
  90 }\r
  91 \r
  92 //______________________________________________________________________________\r
  93 AliAnaOmegaToPi0Gamma & AliAnaOmegaToPi0Gamma::operator = (const AliAnaOmegaToPi0Gamma & ex)\r
  94 {\r
  95  // assignment operator\r
  96 \r
  97  if(this == &ex)return *this;\r
  98    ((AliAnaPartCorrBaseClass *)this)->operator=(ex);\r
  99    fInputAODGamma = new TClonesArray(*ex.fInputAODGamma);\r
 100    fInputAODPi0 = new TClonesArray(*ex.fInputAODPi0);\r
 101 //   fInputAODGamma=ex.fInputAODGamma;\r
 102 //   fInputAODPi0=ex.fInputAODPi0;\r
 103    fInputAODGammaName = ex.fInputAODGammaName;\r
 104    fEventsList = ex.fEventsList;\r
 105 \r
 106    fVtxZCut=ex.fVtxZCut;\r
 107    fCent=ex.fCent;\r
 108    fRp=ex.fRp;\r
 109    fBadChDist=ex.fBadChDist;\r
 110    fNVtxZBin=ex.fNVtxZBin; \r
 111    fNCentBin=ex.fNCentBin;\r
 112    fNRpBin=ex.fNRpBin;\r
 113    fNBadChDistBin=ex.fNBadChDistBin;\r
 114    fNpid=ex.fNpid;\r
 115 \r
 116    fNmaxMixEv =ex.fNmaxMixEv;\r
 117    fPi0Mass=ex.fPi0Mass;\r
 118    fPi0MassWindow=ex.fPi0MassWindow;\r
 119 \r
 120    fNbinsPt=ex.fNbinsPt;fPtBegin=ex.fPtBegin;fPtEnd=ex.fPtEnd;\r
 121    fNbinsM=ex.fNbinsM;fMinM=ex.fMinM;fMaxM=ex.fMaxM;\r
 122    \r
 123    fhEtalon=ex.fhEtalon;\r
 124    fRealOmega=ex.fRealOmega;\r
 125    fMixAOmega=ex.fMixAOmega;\r
 126    fMixBOmega=ex.fMixBOmega;\r
 127    fMixCOmega=ex.fMixCOmega;\r
 128    fRealOmega1=ex.fRealOmega1;\r
 129    fMixAOmega1=ex.fMixAOmega1;\r
 130    fMixBOmega1=ex.fMixBOmega1;\r
 131    fMixCOmega1=ex.fMixCOmega1;\r
 132    fRealOmega2=ex.fRealOmega2;\r
 133    fMixAOmega2=ex.fMixAOmega2;\r
 134    fMixBOmega2=ex.fMixBOmega2;\r
 135    fMixCOmega2=ex.fMixCOmega2;\r
 136    fhOmegaPriPt=ex.fhOmegaPriPt;\r
 137   return *this;\r
 138         \r
 139 }\r
 140 \r
 141 //______________________________________________________________________________\r
 142 AliAnaOmegaToPi0Gamma::~AliAnaOmegaToPi0Gamma() {\r
 143 \r
 144   //dtor\r
 145   if(fInputAODGamma){\r
 146     fInputAODGamma->Clear();\r
 147     delete fInputAODGamma;\r
 148   }\r
 149   \r
 150   if(fInputAODPi0){\r
 151     fInputAODPi0->Clear();\r
 152     delete fInputAODPi0;\r
 153   }  \r
 154 \r
 155   if(fEventsList){\r
 156      for(Int_t i=0;i<fNVtxZBin;i++){\r
 157         for(Int_t j=0;j<fNCentBin;j++){\r
 158            for(Int_t k=0;k<fNRpBin;k++){\r
 159                fEventsList[i*fNCentBin*fNRpBin+j*fNRpBin+k]->Clear();\r
 160                delete fEventsList[i*fNCentBin*fNRpBin+j*fNRpBin+k];\r
 161            }\r
 162         }\r
 163      }\r
 164   }\r
 165   delete [] fEventsList;\r
 166   fEventsList=0;\r
 167 }\r
 168 \r
 169 //______________________________________________________________________________\r
 170 void AliAnaOmegaToPi0Gamma::InitParameters()\r
 171 {\r
 172 //Init parameters when first called the analysis\r
 173 //Set default parameters\r
 174   fInputAODGammaName = "PhotonsPHOS";\r
 175   fNmaxMixEv = 4;\r
 176 \r
 177   fPi0Mass=0.13498;\r
 178   fPi0MassWindow=0.015;\r
 179   fNbinsPt=200;\r
 180   fPtBegin=0;\r
 181   fPtEnd=20;\r
 182   fNbinsM=100;  //mass bin number, min and max\r
 183   fMinM=0;\r
 184   fMaxM=1;\r
 185 \r
 186   fNVtxZBin=1;\r
 187   fNCentBin=1;\r
 188   fNRpBin=1;\r
 189   fNBadChDistBin=1;\r
 190   fNpid=3;\r
 191 \r
 192   fVtxZCut = new Double_t [fNVtxZBin];\r
 193   fVtxZCut[0]=30;\r
 194   //fVtxZCut[0]=10 ; fVtxZCut[1]=20 ; fVtxZCut[2]=30 ;// fVtxZCut[3]=40 ;;\r
 195 \r
 196   fCent=new Double_t[fNCentBin];\r
 197   fCent[0]=0;\r
 198 \r
 199   fRp=new Double_t[fNRpBin];\r
 200   fRp[0]=0;\r
 201 \r
 202   fBadChDist=new Int_t [fNBadChDistBin];\r
 203 //  fBadChDist[0]=0; fBadChDist[1]=1; fBadChDist[2]=2;\r
 204    fBadChDist[0]=0; \r
 205 \r
 206   fhEtalon = new TH2F("hEtalon","Histo with binning parameters", fNbinsPt,fPtBegin, fPtEnd,fNbinsM,fMinM,fMaxM) ;\r
 207   fhEtalon->SetXTitle("P_{T} (GeV)") ;\r
 208   fhEtalon->SetYTitle("m_{inv} (GeV)") ;\r
 209 \r
 210 }\r
 211 \r
 212 \r
 213 //______________________________________________________________________________\r
 214 TList * AliAnaOmegaToPi0Gamma::GetCreateOutputObjects()\r
 215 {  \r
 216  // store them in fOutputContainer\r
 217  fEventsList = new TList*[fNVtxZBin*fNCentBin*fNRpBin];\r
 218  for(Int_t i=0;i<fNVtxZBin;i++){\r
 219      for(Int_t j=0;j<fNCentBin;j++){\r
 220          for(Int_t k=0;k<fNRpBin;k++){\r
 221              fEventsList[i*fNCentBin*fNRpBin+j*fNRpBin+k]=new TList();\r
 222          }\r
 223      }\r
 224  }\r
 225 \r
 226         \r
 227  TList * outputContainer = new TList() ; \r
 228  outputContainer->SetName(GetName());\r
 229  char key[255] ;\r
 230  char title[255] ;\r
 231  const char * detector= fInputAODGammaName.Data();\r
 232  Int_t ndim=fNVtxZBin*fNCentBin*fNRpBin*fNBadChDistBin*fNpid;\r
 233 \r
 234  fRealOmega =new TH2F*[ndim];\r
 235  fMixAOmega =new TH2F*[ndim];\r
 236  fMixBOmega =new TH2F*[ndim];\r
 237  fMixCOmega =new TH2F*[ndim];\r
 238  \r
 239  fRealOmega1 =new TH2F*[ndim];\r
 240  fMixAOmega1 =new TH2F*[ndim];\r
 241  fMixBOmega1 =new TH2F*[ndim];\r
 242  fMixCOmega1 =new TH2F*[ndim];\r
 243 \r
 244  fRealOmega2 =new TH2F*[ndim];\r
 245  fMixAOmega2 =new TH2F*[ndim];\r
 246  fMixBOmega2 =new TH2F*[ndim];\r
 247  fMixCOmega2 =new TH2F*[ndim];\r
 248 \r
 249  for(Int_t i=0;i<fNVtxZBin;i++){\r
 250      for(Int_t j=0;j<fNCentBin;j++){\r
 251          for(Int_t k=0;k<fNRpBin;k++){ //at event level\r
 252              Int_t idim=i*fNCentBin*fNRpBin+j*fNRpBin+k;\r
 253              for(Int_t ipid=0;ipid<fNpid;ipid++){\r
 254                  for(Int_t idist=0;idist<fNBadChDistBin;idist++){ //at particle level\r
 255 \r
 256                      Int_t index=idim*fNpid*fNBadChDistBin+ipid*fNBadChDistBin+idist;\r
 257 \r
 258                      sprintf(key,"RPi0Gamma_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 259                      sprintf(title, "%s Real Pi0GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 260                      fhEtalon->Clone(key);\r
 261                      fRealOmega[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 262                      fRealOmega[index]->SetName(key) ;\r
 263                      fRealOmega[index]->SetTitle(title);\r
 264                      outputContainer->Add(fRealOmega[index]);\r
 265 \r
 266                      sprintf(key,"MAPi0Gamma_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 267                      sprintf(title, "%s MixA Pi0GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 268                      fhEtalon->Clone(key);\r
 269                      fMixAOmega[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 270                      fMixAOmega[index]->SetName(key) ;\r
 271                      fMixAOmega[index]->SetTitle(title);\r
 272                      outputContainer->Add(fMixAOmega[index]);\r
 273 \r
 274                      sprintf(key,"MBPi0Gamma_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 275                      sprintf(title, "%s MixB Pi0GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 276                      fhEtalon->Clone(key);\r
 277                      fMixBOmega[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 278                      fMixBOmega[index]->SetName(key) ;\r
 279                      fMixBOmega[index]->SetTitle(title);\r
 280                      outputContainer->Add(fMixBOmega[index]);\r
 281 \r
 282                      sprintf(key,"MCPi0Gamma_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 283                      sprintf(title, "%s MixC Pi0GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 284                      fhEtalon->Clone(key);\r
 285                      fMixCOmega[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 286                      fMixCOmega[index]->SetName(key) ;\r
 287                      fMixCOmega[index]->SetTitle(title);\r
 288                      outputContainer->Add(fMixCOmega[index]);\r
 289 \r
 290                      sprintf(key,"RPi0Gamma1_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 291                      sprintf(title, "%s Real Pi0(A<0.7)GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 292                      fhEtalon->Clone(key);\r
 293                      fRealOmega1[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 294                      fRealOmega1[index]->SetName(key) ;\r
 295                      fRealOmega1[index]->SetTitle(title);\r
 296                      outputContainer->Add(fRealOmega1[index]);\r
 297 \r
 298                      sprintf(key,"MAPi0Gamma1_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 299                      sprintf(title, "%s MixA Pi0(A<0.7)GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 300                      fhEtalon->Clone(key);\r
 301                      fMixAOmega1[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 302                      fMixAOmega1[index]->SetName(key) ;\r
 303                      fMixAOmega1[index]->SetTitle(title);\r
 304                      outputContainer->Add(fMixAOmega1[index]);\r
 305 \r
 306                      sprintf(key,"MBPi0Gamma1_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 307                      sprintf(title, "%s MixB Pi0(A<0.7)GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 308                      fhEtalon->Clone(key);\r
 309                      fMixBOmega1[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 310                      fMixBOmega1[index]->SetName(key) ;\r
 311                      fMixBOmega1[index]->SetTitle(title);\r
 312                      outputContainer->Add(fMixBOmega1[index]);\r
 313 \r
 314                      sprintf(key,"MCPi0Gamma1_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 315                      sprintf(title, "%s MixC Pi0(A<0.7)GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 316                      fhEtalon->Clone(key);\r
 317                      fMixCOmega1[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 318                      fMixCOmega1[index]->SetName(key) ;\r
 319                      fMixCOmega1[index]->SetTitle(title);\r
 320                      outputContainer->Add(fMixCOmega1[index]);\r
 321 \r
 322                      sprintf(key,"RPi0Gamma2_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 323                      sprintf(title, "%s Real Pi0(A<0.8)GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 324                      fhEtalon->Clone(key);\r
 325                      fRealOmega2[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 326                      fRealOmega2[index]->SetName(key) ;\r
 327                      fRealOmega2[index]->SetTitle(title);\r
 328                      outputContainer->Add(fRealOmega2[index]);\r
 329 \r
 330                      sprintf(key,"MAPi0Gamma2_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 331                      sprintf(title, "%s MixA Pi0(A<0.8)GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 332                      fhEtalon->Clone(key);\r
 333                      fMixAOmega2[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 334                      fMixAOmega2[index]->SetName(key) ;\r
 335                      fMixAOmega2[index]->SetTitle(title);\r
 336                      outputContainer->Add(fMixAOmega2[index]);\r
 337 \r
 338                      sprintf(key,"MBPi0Gamma2_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 339                      sprintf(title, "%s MixB Pi0(A<0.8)GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 340                      fhEtalon->Clone(key);\r
 341                      fMixBOmega2[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 342                      fMixBOmega2[index]->SetName(key) ;\r
 343                      fMixBOmega2[index]->SetTitle(title);\r
 344                      outputContainer->Add(fMixBOmega2[index]);\r
 345 \r
 346                      sprintf(key,"MCPi0Gamma2_Vz%dC%dRp%dPid%dDist%d",i,j,k,ipid,idist);\r
 347                      sprintf(title, "%s MixC Pi0(A<0.8)GammaIVM vz_%2.1f_ct_%2.1f_Rp_%2.1f_pid_%d_dist_%d",detector,fVtxZCut[i],fCent[j],fRp[k],ipid,fBadChDist[idist]);\r
 348                      fhEtalon->Clone(key);\r
 349                      fMixCOmega2[index]=(TH2F*)fhEtalon->Clone(key) ;\r
 350                      fMixCOmega2[index]->SetName(key) ;\r
 351                      fMixCOmega2[index]->SetTitle(title);\r
 352                      outputContainer->Add(fMixCOmega2[index]);\r
 353                  }\r
 354              }\r
 355          }\r
 356      }  \r
 357  }\r
 358 \r
 359  if(IsDataMC()){\r
 360     sprintf(key, "%sOmegaPri",detector);\r
 361     sprintf(title,"primary #omega in %s",detector);\r
 362     fhOmegaPriPt=new TH1F(key, title,fNbinsPt,fPtBegin,fPtEnd);\r
 363     fhOmegaPriPt->GetXaxis()->SetTitle("P_{T}");\r
 364     fhOmegaPriPt->GetYaxis()->SetTitle("dN/P_{T}");\r
 365     outputContainer->Add(fhOmegaPriPt);\r
 366  }\r
 367 \r
 368  delete fhEtalon;\r
 369  return outputContainer;\r
 370 }\r
 371 \r
 372 //______________________________________________________________________________\r
 373 void AliAnaOmegaToPi0Gamma::Print(const Option_t * /*opt*/) const\r
 374 {\r
 375   //Print some relevant parameters set in the analysis\r
 376   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;\r
 377   AliAnaPartCorrBaseClass::Print(" ");\r
 378   printf("Omega->pi0+gamma->3gamma\n");\r
 379   printf("Cuts at event level:            \n");\r
 380   printf("Bins of vertex Z:                     %d \n", fNVtxZBin);\r
 381   printf("Bins of centrality:                   %d \n",fNCentBin);\r
 382   printf("Bins of Reaction plane:               %d\n",fNRpBin);\r
 383   printf("Cuts at AOD particle level:\n");\r
 384   printf("Number of PID:                        %d \n", fNpid);\r
 385   printf("Number of DistToBadChannel cuts:      %d\n", fNBadChDistBin);\r
 386   printf("number of events buffer to be mixed:  %d\n",fNmaxMixEv);\r
 387 } \r
 388 \r
 389 //______________________________________________________________________________\r
 390 void AliAnaOmegaToPi0Gamma::MakeAnalysisFillHistograms() \r
 391 {\r
 392  //fill the MC AOD if needed first\r
 393  //-----------\r
 394  //need to be further implemented\r
 395  AliStack * stack = 0x0;\r
 396 // TParticle * primary = 0x0;\r
 397  TClonesArray * mcparticles0 = 0x0;\r
 398  TClonesArray * mcparticles1 = 0x0;\r
 399  AliAODMCParticle * aodprimary = 0x0;\r
 400  Int_t pdg=0;\r
 401  Double_t pt=0;\r
 402  Double_t eta=0;\r
 403 \r
 404  if(IsDataMC()){\r
 405      if(GetReader()->ReadStack()){\r
 406         stack =  GetMCStack() ;\r
 407         if(!stack)\r
 408         printf("AliAnaAcceptance::MakeAnalysisFillHistograms() - There is no stack!\n");\r
 409         for(Int_t i=0 ; i<stack->GetNtrack(); i++){\r
 410             TParticle * prim = stack->Particle(i) ;\r
 411             pdg = prim->GetPdgCode() ;\r
 412             eta=prim->Eta();\r
 413             pt=prim->Pt();\r
 414             if(TMath::Abs(eta)<0.5) {\r
 415                if(pdg==223) fhOmegaPriPt->Fill(pt);\r
 416             }\r
 417         }\r
 418       }\r
 419       else if(GetReader()->ReadAODMCParticles()){\r
 420           //Get the list of MC particles\r
 421           mcparticles0 = GetReader()->GetAODMCParticles(0);\r
 422           if(!mcparticles0 && GetDebug() > 0)     {\r
 423              printf("AliAnaAcceptance::MakeAnalysisFillHistograms() -  Standard MCParticles not available!\n");\r
 424            }\r
 425            if(GetReader()->GetSecondInputAODTree()){\r
 426                mcparticles1 = GetReader()->GetAODMCParticles(1);\r
 427                if(!mcparticles1 && GetDebug() > 0)     {\r
 428                    printf("AliAnaAcceptance::MakeAnalysisFillHistograms() -  Second input MCParticles not available!\n");\r
 429                 }\r
 430            }\r
 431            for(Int_t i=0;i<mcparticles0->GetEntries();i++){\r
 432                aodprimary =(AliAODMCParticle*)mcparticles0->At(i);\r
 433                pdg = aodprimary->GetPdgCode() ;\r
 434                eta=aodprimary->Eta();\r
 435                pt=aodprimary->Pt();\r
 436                if(TMath::Abs(eta)<0.5) {\r
 437                   if(pdg==223) fhOmegaPriPt->Fill(pt);\r
 438                 }\r
 439 \r
 440            }\r
 441      }\r
 442  }// is data and MC\r
 443   \r
 444 \r
 445  //process event from AOD brach \r
 446  //extract pi0, eta and omega analysis\r
 447  Int_t iRun=(GetReader()->GetInputEvent())->GetRunNumber() ;\r
 448  if(IsBadRun(iRun)) return ;    \r
 449 \r
 450  //vertex z\r
 451  Double_t vert[]={0,0,0} ;\r
 452  GetReader()->GetVertex(vert);\r
 453  Int_t curEventBin =0;\r
 454  Int_t ivtxzbin = -1;\r
 455  for(Int_t i=0;i<fNVtxZBin;i++){\r
 456  //   if(TMath::Abs(vert[2])<fVtxZCut[i] && TMath::Abs(vert[2])>=(fVtxZCut[i]-10)) ivtxzbin=i;\r
 457      if(TMath::Abs(vert[2])<fVtxZCut[i]) ivtxzbin=i;\r
 458  }\r
 459  //centrality\r
 460  Int_t icentbin=0;\r
 461 \r
 462  //reaction plane\r
 463  Int_t irpbin=0;\r
 464 \r
 465  if(ivtxzbin==-1) return; \r
 466  curEventBin = ivtxzbin*fNCentBin*fNRpBin + icentbin*fNRpBin + irpbin;\r
 467 \r
 468   fInputAODGamma = GetAODBranch(fInputAODGammaName); //photon array\r
 469  //fInputAODGamma = (TClonesArray *) GetReader()->GetOutputEvent()->FindListObject(fInputAODGammaName); //photon array\r
 470  Int_t nphotons = fInputAODGamma->GetEntries();\r
 471 \r
 472  fInputAODPi0 = (TClonesArray*)GetInputAODBranch();  //pi0 array\r
 473  Int_t npi0s = fInputAODPi0 ->GetEntries();\r
 474  \r
 475  if(nphotons<3)return; //for pi0, eta and omega->pi0+gamma->3gamma reconstruction\r
 476 \r
 477  //if((nphotons<3) || (npi0s<=0)) return;\r
 478  //reconstruction of omega(782)->pi0+gamma->3gamma\r
 479  //loop for pi0 and photon\r
 480  if(GetDebug() > 0) printf("omega->pi0+gamma->3gamma invariant mass analysis ! \n");\r
 481  for(Int_t i=0;i<npi0s;i++){\r
 482      AliAODPWG4Particle * pi0 = (AliAODPWG4Particle*) (fInputAODPi0->At(i)) ; //pi0\r
 483      TLorentzVector vpi0(pi0->Px(),pi0->Py(),pi0->Pz(),pi0->E());\r
 484      Int_t lab1=pi0->GetCaloLabel(0);  // photon1 from pi0 decay\r
 485      Int_t lab2=pi0->GetCaloLabel(1);  // photon2 from pi0 decay\r
 486      //for omega->pi0+gamma, it needs at least three photons per event\r
 487 \r
 488      //Get the two decay photons from pi0\r
 489      AliAODPWG4Particle * photon1 =0;\r
 490      AliAODPWG4Particle * photon2 =0;\r
 491      for(Int_t d1=0;d1<nphotons;d1++){\r
 492          for(Int_t d2=0;d2<nphotons;d2++){\r
 493              AliAODPWG4Particle * dp1 = (AliAODPWG4Particle*) (fInputAODGamma->At(d1));\r
 494              AliAODPWG4Particle * dp2 = (AliAODPWG4Particle*) (fInputAODGamma->At(d2));\r
 495              Int_t dlab1=dp1->GetCaloLabel(0);\r
 496              Int_t dlab2=dp2->GetCaloLabel(0);\r
 497              if(dlab1==lab1 && dlab2==lab2){\r
 498                 photon1=dp1;\r
 499                 photon2=dp2;\r
 500              }\r
 501              else continue;\r
 502          }\r
 503     }\r
 504 \r
 505     //caculate the asy and dist of the two photon from pi0 decay\r
 506     TLorentzVector dph1(photon1->Px(),photon1->Py(),photon1->Pz(),photon1->E());\r
 507     TLorentzVector dph2(photon2->Px(),photon2->Py(),photon2->Pz(),photon2->E());\r
 508 \r
 509     Double_t pi0asy= TMath::Abs(dph1.E()-dph2.E())/(dph1.E()+dph2.E());\r
 510 //    Double_t phi1=dph1.Phi();\r
 511 //    Double_t phi2=dph2.Phi();\r
 512 //    Double_t eta1=dph1.Eta();\r
 513 //    Double_t eta2=dph2.Eta();\r
 514 //    Double_t pi0dist=TMath::Sqrt((phi1-phi2)*(phi1-phi2)+(eta1-eta2)*(eta1-eta2));\r
 515 \r
 516 \r
 517     if(pi0->GetPdg()==111  && nphotons>2 && npi0s\r
 518            && TMath::Abs(vpi0.M()-fPi0Mass)<fPi0MassWindow) { //pi0 candidates\r
 519 \r
 520          //avoid the double counting\r
 521          Int_t * dc1= new Int_t[nphotons];\r
 522          Int_t * dc2= new Int_t[nphotons];\r
 523          Int_t index1=0;\r
 524          Int_t index2=0;\r
 525          for(Int_t k=0;k<i;k++){\r
 526              AliAODPWG4Particle * p3=(AliAODPWG4Particle*)(fInputAODPi0->At(k));\r
 527              Int_t lab4=p3->GetCaloLabel(0);\r
 528              Int_t lab5=p3->GetCaloLabel(1);\r
 529              if(lab1==lab4){ dc1[index1]=lab5;  index1++;  }\r
 530              if(lab2==lab5){ dc2[index2]=lab4;  index2++;  }\r
 531          }\r
 532 \r
 533 \r
 534          //loop the pi0 with third gamma\r
 535          for(Int_t j=0;j<nphotons;j++){\r
 536              AliAODPWG4Particle *photon3 = (AliAODPWG4Particle*) (fInputAODGamma->At(j));\r
 537              TLorentzVector dph3(photon3->Px(),photon3->Py(),photon3->Pz(),photon3->E());\r
 538              Int_t lab3=photon3->GetCaloLabel(0);\r
 539 \r
 540              for(Int_t l=0;l<index1;l++) if(lab3==dc1[l]) lab3=-1;\r
 541              for(Int_t l=0;l<index2;l++) if(lab3==dc2[l]) lab3=-1;\r
 542 \r
 543              if(lab3>0 && lab3!=lab1 && lab3!=lab2){\r
 544                 for(Int_t ipid=0;ipid<fNpid;ipid++){\r
 545                     for(Int_t idist=0;idist<fNBadChDistBin;idist++){\r
 546                         Int_t index=curEventBin*fNpid*fNBadChDistBin+ipid*fNBadChDistBin+idist;\r
 547                         if(photon1->IsPIDOK(ipid,AliCaloPID::kPhoton) &&\r
 548                            photon2->IsPIDOK(ipid,AliCaloPID::kPhoton) && \r
 549                            photon3->IsPIDOK(ipid,AliCaloPID::kPhoton) &&\r
 550                            photon1->DistToBad()>=fBadChDist[idist] &&\r
 551                            photon2->DistToBad()>=fBadChDist[idist] &&\r
 552                            photon3->DistToBad()>=fBadChDist[idist] ){\r
 553                          \r
 554                            Double_t pi0gammapt=(vpi0+dph3).Pt();\r
 555                            Double_t pi0gammamass=(vpi0+dph3).M();\r
 556                            //fill the histograms\r
 557                            fRealOmega[index]->Fill(pi0gammapt,pi0gammamass); \r
 558                            if(pi0asy<0.7) fRealOmega1[index]->Fill(pi0gammapt,pi0gammamass);\r
 559                            if(pi0asy<0.8) fRealOmega2[index]->Fill(pi0gammapt,pi0gammamass);\r
 560 //printf("real omega  %d  %2.2f \n", index, pi0gammamass);\r
 561                         }\r
 562                      }\r
 563                 }\r
 564              }\r
 565          }      \r
 566          delete []dc1;\r
 567          delete []dc2;\r
 568 \r
 569          //-------------------------\r
 570          //background analysis\r
 571          //three background\r
 572          // --A   (r1_event1+r2_event1)+r3_event2\r
 573 \r
 574          Int_t nMixed = fEventsList[curEventBin]->GetSize();\r
 575          for(Int_t im=0;im<nMixed;im++){\r
 576              TClonesArray* ev2= (TClonesArray*) (fEventsList[curEventBin]->At(im));\r
 577              for(Int_t mix1=0;mix1<ev2->GetEntries();mix1++){\r
 578                 AliAODPWG4Particle *mix1ph = (AliAODPWG4Particle*) (ev2->At(mix1));     \r
 579                 TLorentzVector vmixph(mix1ph->Px(),mix1ph->Py(),mix1ph->Pz(),mix1ph->E());\r
 580                 for(Int_t ipid=0;ipid<fNpid;ipid++){\r
 581                    for(Int_t idist=0;idist<fNBadChDistBin;idist++){\r
 582                        Int_t index=curEventBin*fNpid*fNBadChDistBin+ipid*fNBadChDistBin+idist;\r
 583                        if(photon1->IsPIDOK(ipid,AliCaloPID::kPhoton)&&\r
 584                           photon2->IsPIDOK(ipid,AliCaloPID::kPhoton)&&\r
 585                           mix1ph->IsPIDOK(ipid,AliCaloPID::kPhoton) &&\r
 586                           photon1->DistToBad()>=fBadChDist[idist] &&\r
 587                           photon2->DistToBad()>=fBadChDist[idist] &&\r
 588                           mix1ph->DistToBad()>=fBadChDist[idist] ){\r
 589 \r
 590                           Double_t pi0gammapt=(vpi0+vmixph).Pt();\r
 591                           Double_t pi0gammamass=(vpi0+vmixph).M();\r
 592                           //fill the histograms\r
 593                           fMixAOmega[index]->Fill(pi0gammapt,pi0gammamass);\r
 594                           if(pi0asy<0.7)fMixAOmega1[index]->Fill(pi0gammapt,pi0gammamass);\r
 595                           if(pi0asy<0.8)fMixAOmega2[index]->Fill(pi0gammapt,pi0gammamass);\r
 596 //printf("mix A  %d  %2.2f \n", index, pi0gammamass);\r
 597                           \r
 598                        }\r
 599                    }\r
 600                 }\r
 601              }\r
 602          }\r
 603      }\r
 604  }\r
 605 \r
 606  //\r
 607  // --B   (r1_event1+r2_event2)+r3_event2\r
 608  //\r
 609  for(Int_t i=0;i<nphotons;i++){\r
 610      AliAODPWG4Particle *ph1 = (AliAODPWG4Particle*) (fInputAODGamma->At(i)); \r
 611      TLorentzVector vph1(ph1->Px(),ph1->Py(),ph1->Pz(),ph1->E());\r
 612 \r
 613      Int_t nMixed = fEventsList[curEventBin]->GetSize();\r
 614      for(Int_t ie=0;ie<nMixed;ie++){\r
 615          TClonesArray* ev2= (TClonesArray*) (fEventsList[curEventBin]->At(ie));\r
 616          for(Int_t mix1=0;mix1<ev2->GetEntries();mix1++){\r
 617              AliAODPWG4Particle *ph2 = (AliAODPWG4Particle*) (ev2->At(mix1));\r
 618              TLorentzVector vph2(ph2->Px(),ph2->Py(),ph2->Pz(),ph2->E());\r
 619              Double_t pi0asy = TMath::Abs(vph1.E()-vph2.E())/(vph1.E()+vph2.E());            \r
 620              Double_t pi0mass=(vph1+vph2).M();\r
 621                 \r
 622              if(TMath::Abs(pi0mass-fPi0Mass)<fPi0MassWindow){//for pi0 selection\r
 623                 for(Int_t mix2=(mix1+1);mix2<ev2->GetEntries();mix2++){\r
 624                     AliAODPWG4Particle *ph3 = (AliAODPWG4Particle*) (ev2->At(mix2));\r
 625                     TLorentzVector vph3(ph3->Px(),ph3->Py(),ph3->Pz(),ph3->E());\r
 626                     \r
 627                     Double_t pi0gammapt=(vph1+vph2+vph3).Pt();\r
 628                     Double_t pi0gammamass=(vph1+vph2+vph3).M(); \r
 629                     for(Int_t ipid=0;ipid<fNpid;ipid++){\r
 630                          for(Int_t idist=0;idist<fNBadChDistBin;idist++){\r
 631                             Int_t index=curEventBin*fNpid*fNBadChDistBin+ipid*fNBadChDistBin+idist;\r
 632                             if(ph1->IsPIDOK(ipid,AliCaloPID::kPhoton) &&\r
 633                                ph2->IsPIDOK(ipid,AliCaloPID::kPhoton) &&\r
 634                                ph3->IsPIDOK(ipid,AliCaloPID::kPhoton) &&\r
 635                                ph1->DistToBad()>=fBadChDist[idist] &&\r
 636                                ph2->DistToBad()>=fBadChDist[idist] &&\r
 637                                ph3->DistToBad()>=fBadChDist[idist] ){\r
 638                                //fill histograms\r
 639                                fMixBOmega[index]->Fill(pi0gammapt,pi0gammamass);\r
 640                                if(pi0asy<0.7) fMixBOmega1[index]->Fill(pi0gammapt,pi0gammamass);\r
 641                                if(pi0asy<0.8) fMixBOmega2[index]->Fill(pi0gammapt,pi0gammamass);\r
 642 //printf("mix B  %d  %2.2f \n", index, pi0gammamass);\r
 643                             }\r
 644                         }                   \r
 645                     }\r
 646                 }\r
 647                 \r
 648                 //\r
 649                 // --C   (r1_event1+r2_event2)+r3_event3\r
 650                 //\r
 651                 for(Int_t je=(ie+1);je<nMixed;je++){\r
 652                     TClonesArray* ev3= (TClonesArray*) (fEventsList[curEventBin]->At(je));\r
 653                     for(Int_t mix3=0;mix3<ev3->GetEntries();mix3++){\r
 654                         AliAODPWG4Particle *ph3 = (AliAODPWG4Particle*) (ev3->At(mix3));\r
 655                         TLorentzVector vph3(ph3->Px(),ph3->Py(),ph3->Pz(),ph3->E());\r
 656                         \r
 657                         Double_t pi0gammapt=(vph1+vph2+vph3).Pt();\r
 658                         Double_t pi0gammamass=(vph1+vph2+vph3).M();\r
 659                         \r
 660                         for(Int_t ipid=0;ipid<fNpid;ipid++){\r
 661                             for(Int_t idist=0;idist<fNBadChDistBin;idist++){\r
 662                                 Int_t index=curEventBin*fNpid*fNBadChDistBin+ipid*fNBadChDistBin+idist;\r
 663                                 if(ph1->IsPIDOK(ipid,AliCaloPID::kPhoton) &&\r
 664                                    ph2->IsPIDOK(ipid,AliCaloPID::kPhoton) &&\r
 665                                    ph3->IsPIDOK(ipid,AliCaloPID::kPhoton) &&\r
 666                                    ph1->DistToBad()>=fBadChDist[idist] &&\r
 667                                    ph2->DistToBad()>=fBadChDist[idist] &&\r
 668                                    ph3->DistToBad()>=fBadChDist[idist] ){\r
 669                                    //fill histograms\r
 670                                   fMixCOmega[index]->Fill(pi0gammapt,pi0gammamass);\r
 671                                   if(pi0asy<0.7) fMixCOmega1[index]->Fill(pi0gammapt,pi0gammamass);\r
 672                                   if(pi0asy<0.8) fMixCOmega2[index]->Fill(pi0gammapt,pi0gammamass);\r
 673 //printf("mix C  %d  %2.2f \n", index, pi0gammamass);\r
 674                                }\r
 675                            }\r
 676                         }\r
 677                     }\r
 678                 }\r
 679              } //for pi0 selecton               \r
 680          }\r
 681      }\r
 682  }\r
 683 \r
 684  //event buffer \r
 685  TClonesArray *currentEvent = new TClonesArray(*fInputAODGamma);\r
 686  if(currentEvent->GetEntriesFast()>0){\r
 687         fEventsList[curEventBin]->AddFirst(currentEvent) ;\r
 688         currentEvent=0 ; \r
 689         if(fEventsList[curEventBin]->GetSize()>=fNmaxMixEv) {\r
 690             TClonesArray * tmp = (TClonesArray*) (fEventsList[curEventBin]->Last()) ;\r
 691             fEventsList[curEventBin]->RemoveLast() ;\r
 692             delete tmp ;\r
 693         }\r
 694  }\r
 695  else{ \r
 696        delete currentEvent ;\r
 697        currentEvent=0 ;\r
 698  }\r
 699  \r
 700 }\r
 701 \r
 702 //______________________________________________________________________________\r
 703 void AliAnaOmegaToPi0Gamma::ReadHistograms(TList * outputList)\r
 704 {\r
 705  //read the histograms \r
 706  //for the finalization of the terminate analysis\r
 707 \r
 708  Int_t index = outputList->IndexOf(outputList->FindObject(GetAddedHistogramsStringToName()+"RPi0Gamma_Vz0C0Rp0Pid0Dist0"));\r
 709 \r
 710   Int_t ndim=fNVtxZBin*fNCentBin*fNRpBin*fNBadChDistBin*fNpid;\r
 711 \r
 712  if(!fRealOmega) fRealOmega =new TH2F*[ndim];\r
 713  if(!fMixAOmega) fMixAOmega =new TH2F*[ndim];\r
 714  if(!fMixBOmega) fMixBOmega =new TH2F*[ndim];\r
 715  if(!fMixCOmega) fMixCOmega =new TH2F*[ndim];\r
 716 \r
 717  if(!fRealOmega1) fRealOmega1 =new TH2F*[ndim];\r
 718  if(!fMixAOmega1) fMixAOmega1 =new TH2F*[ndim];\r
 719  if(!fMixBOmega1) fMixBOmega1 =new TH2F*[ndim];\r
 720  if(!fMixCOmega1) fMixCOmega1 =new TH2F*[ndim];\r
 721 \r
 722  if(!fRealOmega2) fRealOmega2 =new TH2F*[ndim];\r
 723  if(!fMixAOmega2) fMixAOmega2 =new TH2F*[ndim];\r
 724  if(!fMixBOmega2) fMixBOmega2 =new TH2F*[ndim];\r
 725  if(!fMixCOmega2) fMixCOmega2 =new TH2F*[ndim];\r
 726 \r
 727   for(Int_t i=0;i<fNVtxZBin;i++){\r
 728      for(Int_t j=0;j<fNCentBin;j++){\r
 729          for(Int_t k=0;k<fNRpBin;k++){ //at event level\r
 730              Int_t idim=i*fNCentBin*fNRpBin+j*fNRpBin+k;\r
 731              for(Int_t ipid=0;ipid<fNpid;ipid++){ \r
 732                 for(Int_t idist=0;idist<fNBadChDistBin;idist++){ //at particle\r
 733                     Int_t ind=idim*fNpid*fNBadChDistBin+ipid*fNBadChDistBin+idist;\r
 734                     fRealOmega[ind]= (TH2F*) outputList->At(index++);\r
 735                     fMixAOmega[ind]= (TH2F*) outputList->At(index++);\r
 736                     fMixBOmega[ind]= (TH2F*) outputList->At(index++);\r
 737                     fMixCOmega[ind]= (TH2F*) outputList->At(index++);\r
 738 \r
 739                     fRealOmega1[ind]= (TH2F*) outputList->At(index++);\r
 740                     fMixAOmega1[ind]= (TH2F*) outputList->At(index++);\r
 741                     fMixBOmega1[ind]= (TH2F*) outputList->At(index++);\r
 742                     fMixCOmega1[ind]= (TH2F*) outputList->At(index++);\r
 743 \r
 744                     fRealOmega2[ind]= (TH2F*) outputList->At(index++);\r
 745                     fMixAOmega2[ind]= (TH2F*) outputList->At(index++);\r
 746                     fMixBOmega2[ind]= (TH2F*) outputList->At(index++);\r
 747                     fMixCOmega2[ind]= (TH2F*) outputList->At(index++);\r
 748                     \r
 749                  \r
 750                 }\r
 751               }\r
 752           }\r
 753       }\r
 754   }\r
 755   \r
 756   if(IsDataMC()){\r
 757      fhOmegaPriPt  = (TH1F*)  outputList->At(index++);\r
 758   }\r
 759 \r
 760 }\r
 761 \r
 762 //______________________________________________________________________________\r
 763 void AliAnaOmegaToPi0Gamma::Terminate(TList * outputList) \r
 764 {\r
 765 // //Do some calculations and plots from the final histograms.\r
 766   if(GetDebug() >= 0) printf("AliAnaOmegaToPi0Gamma::Terminate() \n");\r
 767   ReadHistograms(outputList);\r
 768   char cvs1[128];  \r
 769   sprintf(cvs1, "Neutral_%s_IVM",fInputAODGammaName.Data());\r
 770 \r
 771   TCanvas * cvsIVM = new TCanvas(cvs1, cvs1, 400, 10, 600, 700) ;\r
 772   cvsIVM->Divide(2, 2);\r
 773 \r
 774   cvsIVM->cd(1);\r
 775   TH2F * h2Real= new TH2F(*fRealOmega[0]);\r
 776   h2Real->GetXaxis()->SetRangeUser(4,6);\r
 777   TH1F * hRealOmega = (TH1F*) h2Real->ProjectionY();\r
 778   hRealOmega->SetTitle("RealPi0Gamma 4<pt<6");\r
 779   hRealOmega->SetLineColor(2);\r
 780   hRealOmega->Draw();\r
 781 \r
 782   cvsIVM->cd(2);\r
 783   TH2F *h2MixA= new TH2F(*fMixAOmega[0]);\r
 784   h2MixA->GetXaxis()->SetRangeUser(4,6);\r
 785   TH1F * hMixAOmega = (TH1F*) h2MixA->ProjectionY();\r
 786   hMixAOmega->SetTitle("MixA 4<pt<6");\r
 787   hMixAOmega->SetLineColor(2);\r
 788   hMixAOmega->Draw();\r
 789 \r
 790   cvsIVM->cd(3);\r
 791   TH2F * h2MixB=new TH2F(*fMixBOmega[0]);\r
 792   h2MixB->GetXaxis()->SetRangeUser(4,6);\r
 793   TH1F * hMixBOmega = (TH1F*) h2MixB->ProjectionY();\r
 794   hMixBOmega->SetTitle("MixB 4<pt<6");\r
 795   hMixBOmega->SetLineColor(2);\r
 796   hMixBOmega->Draw();\r
 797 \r
 798   cvsIVM->cd(4);\r
 799   TH2F *h2MixC= new TH2F(*fMixCOmega[0]);\r
 800   h2MixC->GetXaxis()->SetRangeUser(4,6);\r
 801   TH1F * hMixCOmega = (TH1F*) h2MixC->ProjectionY();\r
 802   hMixCOmega->SetTitle("MixC 4<pt<6");\r
 803   hMixCOmega->SetLineColor(2);\r
 804   hMixCOmega->Draw();\r
 805 \r
 806   char eps[128 ];\r
 807   sprintf(eps,"CVS_%s_IVM.eps",fInputAODGammaName.Data());\r
 808   cvsIVM->Print(eps);\r
 809   cvsIVM->Modified();\r
 810  \r
 811 }\r