]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/PHOSTasks/PHOS_embedding/AliAnalysisTaskPi0DiffEfficiency.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
AliAODEvent::GetHeader() returns AliVHeader
[u/mrichter/AliRoot.git] / PWGGA / PHOSTasks / PHOS_embedding / AliAnalysisTaskPi0DiffEfficiency.cxx
CommitLineData
02f0b8e7 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"
7063ee4d 14#include "THashList.h"
02f0b8e7 15
16#include "AliAODMCParticle.h"
17#include "AliAnalysisManager.h"
18#include "AliMCEventHandler.h"
19#include "AliMCEvent.h"
20#include "AliStack.h"
21#include "AliAnalysisTaskSE.h"
22#include "AliAnalysisTaskPi0DiffEfficiency.h"
23#include "AliCaloPhoton.h"
24#include "AliPHOSGeometry.h"
6e4947dd 25#include "AliPHOSAodCluster.h"
02f0b8e7 26#include "AliPHOSCalibData.h"
27#include "AliAODEvent.h"
28#include "AliAODCaloCluster.h"
29#include "AliAODVertex.h"
30#include "AliESDtrackCuts.h"
31#include "AliLog.h"
32#include "AliPID.h"
33#include "AliCDBManager.h"
34#include "AliCentrality.h"
35
6e4947dd 36// Analysis task to calculate pi0 registration efficiency
37// as a difference between spectra before and after embedding
38// Authors: Dmitri Peressounko
39// Date : Aug.2011
02f0b8e7 40
7063ee4d 41Double_t Scale(Double_t x){
42
6631c482 43// return 1./1.008/1.015*(1.+0.017/(1.+x*x/2./2.)+0.03/(1.+x*x/0.6/0.6)) ;
44 return 1./1.015/1.015*(1.+0.017/(1.+x*x/2./2.)+0.03/(1.+x*x/0.6/0.6)) ;
7063ee4d 45
46
47}
48
02f0b8e7 49ClassImp(AliAnalysisTaskPi0DiffEfficiency)
fc818b37 50//________________________________________________________________________
51AliAnalysisTaskPi0DiffEfficiency::AliAnalysisTaskPi0DiffEfficiency()
52: AliAnalysisTaskSE(),
53 fStack(0),
54 fOutputContainer(0),
55 fPHOSEvent(0),
56 fPHOSEvent1(0),
57 fPHOSEvent2(0),
58 fPHOSCalibData(0),
59 fNonLinCorr(0),
60 fRPfull(0),
61 fRPA(0),
62 fRPC(0),
63 fRPFar(0),
64 fRPAFar(0),
65 fRPCFar(0),
66 fCentrality(0),
67 fCenBin(0),
68 fPHOSGeo(0),
69 fEventCounter(0)
70{
71 // Constructor //should not be used
72 for(Int_t i=0;i<1;i++){
73 for(Int_t j=0;j<10;j++)
74 for(Int_t k=0;k<11;k++)
75 fPHOSEvents[i][j][k]=0 ;
76 }
7bf900cb 77 for(Int_t mod=0; mod<6; mod++) fPHOSBadMap[mod]=0x0 ;
fc818b37 78}
02f0b8e7 79
80//________________________________________________________________________
81AliAnalysisTaskPi0DiffEfficiency::AliAnalysisTaskPi0DiffEfficiency(const char *name)
7063ee4d 82: AliAnalysisTaskSE(name),
83 fStack(0),
84 fOutputContainer(0),
85 fPHOSEvent(0),
02f0b8e7 86 fPHOSEvent1(0),
7063ee4d 87 fPHOSEvent2(0),
88 fPHOSCalibData(0),
89 fNonLinCorr(0),
90 fRPfull(0),
91 fRPA(0),
92 fRPC(0),
93 fRPFar(0),
94 fRPAFar(0),
95 fRPCFar(0),
96 fCentrality(0),
97 fCenBin(0),
98 fPHOSGeo(0),
fc818b37 99 fEventCounter(0)
02f0b8e7 100{
7063ee4d 101 // Constructor
102 for(Int_t i=0;i<1;i++){
103 for(Int_t j=0;j<10;j++)
104 for(Int_t k=0;k<11;k++)
105 fPHOSEvents[i][j][k]=0 ;
106 }
107
108 // Output slots #0 write into a TH1 container
109 DefineOutput(1,TList::Class());
110
111 // Set bad channel map
112 char key[55] ;
113 for(Int_t i=0; i<6; i++){
114 snprintf(key,55,"PHOS_BadMap_mod%d",i) ;
115 fPHOSBadMap[i]=new TH2I(key,"Bad Modules map",64,0.,64.,56,0.,56.) ;
116 }
117 // Initialize the PHOS geometry
118 fPHOSGeo = AliPHOSGeometry::GetInstance("IHEP") ;
119
120 fPHOSCalibData = new AliPHOSCalibData();
121 for(Int_t module=1; module<=5; module++) {
122 for(Int_t column=1; column<=56; column++) {
123 for(Int_t row=1; row<=64; row++) {
124 fPHOSCalibData->SetADCchannelEmc(module,column,row,1.);
125 }
126 }
127 }
02f0b8e7 128
02f0b8e7 129
130
131}
fc818b37 132//________________________________________________________________________
133AliAnalysisTaskPi0DiffEfficiency::AliAnalysisTaskPi0DiffEfficiency(const AliAnalysisTaskPi0DiffEfficiency& a):
134AliAnalysisTaskSE(a),
135 fStack(0),
136 fOutputContainer(0),
137 fPHOSEvent(0),
138 fPHOSEvent1(0),
139 fPHOSEvent2(0),
140 fPHOSCalibData(0),
141 fNonLinCorr(0),
142 fRPfull(0),
143 fRPA(0),
144 fRPC(0),
145 fRPFar(0),
146 fRPAFar(0),
147 fRPCFar(0),
148 fCentrality(0),
149 fCenBin(0),
150 fPHOSGeo(0),
151 fEventCounter(0)
152{ // not implemented
7bf900cb 153 for(Int_t mod=0; mod<6; mod++) fPHOSBadMap[mod]=0x0 ;
fc818b37 154}
02f0b8e7 155//________________________________________________________________________
156void AliAnalysisTaskPi0DiffEfficiency::UserCreateOutputObjects()
157{
158 // Create histograms
159 // Called once
160
161 // ESD histograms
162 if(fOutputContainer != NULL){
163 delete fOutputContainer;
164 }
7063ee4d 165 fOutputContainer = new THashList();
02f0b8e7 166 fOutputContainer->SetOwner(kTRUE);
167
168 //Event selection
169 fOutputContainer->Add(new TH1F("hSelEvents","Event celection", 10,0.,10.)) ;
170
171 //vertex distribution
172 fOutputContainer->Add(new TH1F("hZvertex","Z vertex position", 50,-25.,25.)) ;
173
174 //Centrality
175 fOutputContainer->Add(new TH1F("hCentrality","Event centrality", 100,0.,100.)) ;
176
177 //QA histograms
178 fOutputContainer->Add(new TH2F("hCluM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
179 fOutputContainer->Add(new TH2F("hCluM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
180 fOutputContainer->Add(new TH2F("hCluM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
181
182 Int_t nM = 500;
183 Double_t mMin = 0.0;
184 Double_t mMax = 1.0;
6631c482 185 Int_t nPt = 250;
02f0b8e7 186 Double_t ptMin = 0;
6631c482 187 Double_t ptMax = 25;
7063ee4d 188
6631c482 189 Int_t nPtf = 23;
190 Double_t xPt[24]={0.6,1.,1.5,2.,2.5,3.,3.5,4.,4.5,5.,5.5,6.,7.,8.,9.,10.,12.,14.,16.,18.,20.,22.,24.,25.} ;
7063ee4d 191 Int_t nPhi=10 ;
192 Double_t xPhi[11] ;
193 for(Int_t i=0;i<=10;i++)xPhi[i]=i*0.1*TMath::Pi() ;
194 Int_t nMm=150 ;
195 Double_t xM[201] ;
196 for(Int_t i=0;i<=200;i++){xM[i]=0.0025*i;}
197
02f0b8e7 198
199 char key[55] ;
200 for(Int_t cent=0; cent<6; cent++){
201 //Single photon
cf150629 202
203 fOutputContainer->Add(new TH1F(Form("hPhotAll_DistBad2_cen%d",cent),"All clusters",nPt,ptMin,ptMax));
204 fOutputContainer->Add(new TH1F(Form("hPhotAll_DistBad4_cen%d",cent),"All clusters",nPt,ptMin,ptMax));
205 fOutputContainer->Add(new TH1F(Form("hPhotAll_DistBad6_cen%d",cent),"All clusters",nPt,ptMin,ptMax));
02f0b8e7 206 snprintf(key,55,"hPhotAll_cen%d",cent) ;
207 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
208 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 209 snprintf(key,55,"hPhotAllcore_cen%d",cent) ;
210 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
211 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
212 snprintf(key,55,"hPhotAllwou_cen%d",cent) ;
213 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
214 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 215 snprintf(key,55,"hPhotCPV_cen%d",cent) ;
216 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
6e4947dd 217 snprintf(key,55,"hPhotCPVcore_cen%d",cent) ;
218 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
219 snprintf(key,55,"hPhotCPV2_cen%d",cent) ;
220 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
02f0b8e7 221 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 222 snprintf(key,55,"hPhotCPV2core_cen%d",cent) ;
223 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
224 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 225 snprintf(key,55,"hPhotDisp_cen%d",cent) ;
226 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
227 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 228 snprintf(key,55,"hPhotDispwou_cen%d",cent) ;
229 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
230 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 231 snprintf(key,55,"hPhotDisp2_cen%d",cent) ;
232 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
233 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
cf150629 234 snprintf(key,55,"hPhotDispcore_cen%d",cent) ;
6e4947dd 235 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
7063ee4d 236 snprintf(key,55,"hPhotDisp2core_cen%d",cent) ;
237 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
6e4947dd 238 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 239 snprintf(key,55,"hPhotBoth_cen%d",cent) ;
240 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
241 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 242 snprintf(key,55,"hPhotBothcore_cen%d",cent) ;
243 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
244 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 245 snprintf(key,55,"hPhotBoth2_cen%d",cent) ;
246 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
247 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
248 snprintf(key,55,"hPhotBoth2core_cen%d",cent) ;
249 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
250 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 251
6e4947dd 252 snprintf(key,55,"hNegPhotAll_cen%d",cent) ;
253 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
254 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
255 snprintf(key,55,"hNegPhotAllcore_cen%d",cent) ;
256 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
257 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 258 snprintf(key,55,"hNegPhotCPV_cen%d",cent) ;
259 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
260 snprintf(key,55,"hNegPhotCPVcore_cen%d",cent) ;
261 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
262 snprintf(key,55,"hNegPhotCPV2_cen%d",cent) ;
263 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
264 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
265 snprintf(key,55,"hNegPhotDisp_cen%d",cent) ;
266 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
267 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
268 snprintf(key,55,"hNegPhotDisp2_cen%d",cent) ;
269 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
270 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
cf150629 271 snprintf(key,55,"hNegPhotDispcore_cen%d",cent) ;
6e4947dd 272 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
7063ee4d 273 snprintf(key,55,"hNegPhotDisp2core_cen%d",cent) ;
274 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
6e4947dd 275 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
276 snprintf(key,55,"hNegPhotBoth_cen%d",cent) ;
277 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
278 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
279 snprintf(key,55,"hNegPhotBothcore_cen%d",cent) ;
280 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
281 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 282 snprintf(key,55,"hNegPhotBoth2_cen%d",cent) ;
283 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
284 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
285 snprintf(key,55,"hNegPhotBoth2core_cen%d",cent) ;
286 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
287 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 288
02f0b8e7 289 snprintf(key,55,"hOldMassPtAll_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() ;
6e4947dd 292 snprintf(key,55,"hOldMassPtAllcore_cen%d",cent) ;
293 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
294 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 295 snprintf(key,55,"hOldMassPtCPV_cen%d",cent) ;
296 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
6e4947dd 297 snprintf(key,55,"hOldMassPtCPVcore_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,"hOldMassPtCPV2_cen%d",cent) ;
300 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
02f0b8e7 301 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
302 snprintf(key,55,"hOldMassPtDisp_cen%d",cent) ;
303 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
304 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 305 snprintf(key,55,"hOldMassPtDisp2_cen%d",cent) ;
306 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
307 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
cf150629 308 snprintf(key,55,"hOldMassPtDispcore_cen%d",cent) ;
6e4947dd 309 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
310 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 311 snprintf(key,55,"hOldMassPtDisp2core_cen%d",cent) ;
312 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
313 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 314 snprintf(key,55,"hOldMassPtBoth_cen%d",cent) ;
315 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
316 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 317 snprintf(key,55,"hOldMassPtBothcore_cen%d",cent) ;
318 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
319 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 320 snprintf(key,55,"hOldMassPtBoth2_cen%d",cent) ;
321 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
322 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
323 snprintf(key,55,"hOldMassPtBoth2core_cen%d",cent) ;
324 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
325 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 326
327 snprintf(key,55,"hNewMassPtAll_cen%d",cent) ;
328 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
329 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 330 snprintf(key,55,"hNewMassPtAllcore_cen%d",cent) ;
331 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
332 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 333 snprintf(key,55,"hNewMassPtCPV_cen%d",cent) ;
334 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
335 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 336 snprintf(key,55,"hNewMassPtCPVcore_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,"hNewMassPtCPV2_cen%d",cent) ;
339 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
340 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
341 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 342 snprintf(key,55,"hNewMassPtDisp_cen%d",cent) ;
343 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
344 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 345 snprintf(key,55,"hNewMassPtDisp2_cen%d",cent) ;
346 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
347 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
cf150629 348 snprintf(key,55,"hNewMassPtDispcore_cen%d",cent) ;
6e4947dd 349 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
350 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 351 snprintf(key,55,"hNewMassPtDisp2core_cen%d",cent) ;
352 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
353 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 354 snprintf(key,55,"hNewMassPtBoth_cen%d",cent) ;
355 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
356 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 357 snprintf(key,55,"hNewMassPtBothcore_cen%d",cent) ;
358 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
359 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 360 snprintf(key,55,"hNewMassPtBoth2_cen%d",cent) ;
361 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
362 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
363 snprintf(key,55,"hNewMassPtBoth2core_cen%d",cent) ;
364 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
365 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 366
367 snprintf(key,55,"hNegMassPtAll_cen%d",cent) ;
368 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
369 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
370 snprintf(key,55,"hNegMassPtAllcore_cen%d",cent) ;
371 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
372 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 373 snprintf(key,55,"hNegMassPtCPV_cen%d",cent) ;
374 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
375 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
376 snprintf(key,55,"hNegMassPtCPVcore_cen%d",cent) ;
377 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
378 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
379 snprintf(key,55,"hNegMassPtCPV2_cen%d",cent) ;
380 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
381 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
382 snprintf(key,55,"hNegMassPtDisp_cen%d",cent) ;
383 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
384 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
385 snprintf(key,55,"hNegMassPtDisp2_cen%d",cent) ;
386 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
387 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
cf150629 388 snprintf(key,55,"hNegMassPtDispcore_cen%d",cent) ;
6e4947dd 389 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
390 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 391 snprintf(key,55,"hNegMassPtDisp2core_cen%d",cent) ;
392 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
393 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 394 snprintf(key,55,"hNegMassPtBoth_cen%d",cent) ;
395 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
396 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
397 snprintf(key,55,"hNegMassPtBothcore_cen%d",cent) ;
398 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
399 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 400 snprintf(key,55,"hNegMassPtBoth2_cen%d",cent) ;
401 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
402 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
403 snprintf(key,55,"hNegMassPtBoth2core_cen%d",cent) ;
404 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
405 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 406
02f0b8e7 407
408 snprintf(key,55,"hMassPtAll_cen%d",cent) ;
409 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
410 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 411 snprintf(key,55,"hMassPtAllwou_cen%d",cent) ;
6e4947dd 412 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
413 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 414 snprintf(key,55,"hMassPtAllcore_cen%d",cent) ;
6e4947dd 415 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
416 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 417 snprintf(key,55,"hMassPtCPV_cen%d",cent) ;
418 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
419 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 420 snprintf(key,55,"hMassPtCPVcore_cen%d",cent) ;
421 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
422 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
423 snprintf(key,55,"hMassPtCPV2_cen%d",cent) ;
424 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
425 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 426 snprintf(key,55,"hMassPtDisp_cen%d",cent) ;
427 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
428 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 429 snprintf(key,55,"hMassPtDispwou_cen%d",cent) ;
430 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
431 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 432 snprintf(key,55,"hMassPtDisp2_cen%d",cent) ;
433 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
434 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
cf150629 435 snprintf(key,55,"hMassPtDispcore_cen%d",cent) ;
6e4947dd 436 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
437 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 438 snprintf(key,55,"hMassPtDisp2core_cen%d",cent) ;
439 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
440 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 441 snprintf(key,55,"hMassPtBoth_cen%d",cent) ;
442 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
443 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 444 snprintf(key,55,"hMassPtBothcore_cen%d",cent) ;
445 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
446 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 447 snprintf(key,55,"hMassPtBoth2_cen%d",cent) ;
6e4947dd 448 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 449 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
450 snprintf(key,55,"hMassPtBoth2core_cen%d",cent) ;
6e4947dd 451 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 452 ((TH2F*)fOutputContainer->Last())->Sumw2() ;
453
02f0b8e7 454
455 //Mixed
456 snprintf(key,55,"hMiMassPtAll_cen%d",cent) ;
457 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
6e4947dd 458 snprintf(key,55,"hMiMassPtAllwou_cen%d",cent) ;
459 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 460 snprintf(key,55,"hMiMassPtAllcore_cen%d",cent) ;
461 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
02f0b8e7 462 snprintf(key,55,"hMiMassPtCPV_cen%d",cent) ;
463 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
6e4947dd 464 snprintf(key,55,"hMiMassPtCPVcore_cen%d",cent) ;
465 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
466 snprintf(key,55,"hMiMassPtCPV2_cen%d",cent) ;
467 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
02f0b8e7 468 snprintf(key,55,"hMiMassPtDisp_cen%d",cent) ;
469 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 470 snprintf(key,55,"hMiMassPtDispwou_cen%d",cent) ;
471 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
6e4947dd 472 snprintf(key,55,"hMiMassPtDisp2_cen%d",cent) ;
473 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
cf150629 474 snprintf(key,55,"hMiMassPtDispcore_cen%d",cent) ;
6e4947dd 475 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 476 snprintf(key,55,"hMiMassPtDisp2core_cen%d",cent) ;
477 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
02f0b8e7 478 snprintf(key,55,"hMiMassPtBoth_cen%d",cent) ;
479 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
6e4947dd 480 snprintf(key,55,"hMiMassPtBothcore_cen%d",cent) ;
481 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 482 snprintf(key,55,"hMiMassPtBoth2_cen%d",cent) ;
6e4947dd 483 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 484 snprintf(key,55,"hMiMassPtBoth2core_cen%d",cent) ;
6e4947dd 485 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 486
02f0b8e7 487 snprintf(key,55,"hMCMassPtAll_cen%d",cent) ;
488 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
6e4947dd 489 snprintf(key,55,"hMCMassPtAllwou_cen%d",cent) ;
490 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 491 snprintf(key,55,"hMCMassPtAllcore_cen%d",cent) ;
492 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
02f0b8e7 493 snprintf(key,55,"hMCMassPtCPV_cen%d",cent) ;
494 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
6e4947dd 495 snprintf(key,55,"hMCMassPtCPVcore_cen%d",cent) ;
496 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
497 snprintf(key,55,"hMCMassPtCPV2_cen%d",cent) ;
498 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 499 snprintf(key,55,"hMCMassPtCPV2core_cen%d",cent) ;
500 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
02f0b8e7 501 snprintf(key,55,"hMCMassPtDisp_cen%d",cent) ;
502 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 503 snprintf(key,55,"hMCMassPtDispwou_cen%d",cent) ;
6e4947dd 504 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
cf150629 505 snprintf(key,55,"hMCMassPtDispcore_cen%d",cent) ;
6e4947dd 506 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 507 snprintf(key,55,"hMCMassPtDisp2_cen%d",cent) ;
6e4947dd 508 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 509 snprintf(key,55,"hMCMassPtDisp2core_cen%d",cent) ;
6e4947dd 510 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 511 snprintf(key,55,"hMCMassPtBoth_cen%d",cent) ;
6e4947dd 512 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 513 snprintf(key,55,"hMCMassPtBoth2_cen%d",cent) ;
6e4947dd 514 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 515 snprintf(key,55,"hMCMassPtBothcore_cen%d",cent) ;
6e4947dd 516 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
7063ee4d 517 snprintf(key,55,"hMCMassPtBoth2core_cen%d",cent) ;
6e4947dd 518 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
02f0b8e7 519
520 //Single photon
521 snprintf(key,55,"hMCPhotAll_cen%d",cent) ;
522 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
523 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 524 snprintf(key,55,"hMCPhotAllwou_cen%d",cent) ;
6e4947dd 525 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
526 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 527 snprintf(key,55,"hMCPhotAllcore_cen%d",cent) ;
6e4947dd 528 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
529 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 530 snprintf(key,55,"hMCPhotCPV_cen%d",cent) ;
531 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
532 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 533 snprintf(key,55,"hMCPhotCPVcore_cen%d",cent) ;
534 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
535 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
536 snprintf(key,55,"hMCPhotCPV2_cen%d",cent) ;
537 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
538 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
02f0b8e7 539 snprintf(key,55,"hMCPhotDisp_cen%d",cent) ;
540 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
541 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 542 snprintf(key,55,"hMCPhotDispwou_cen%d",cent) ;
543 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
544 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 545 snprintf(key,55,"hMCPhotDisp2_cen%d",cent) ;
02f0b8e7 546 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
547 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
cf150629 548 snprintf(key,55,"hMCPhotDispcore_cen%d",cent) ;
6e4947dd 549 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
550 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 551 snprintf(key,55,"hMCPhotDisp2core_cen%d",cent) ;
552 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
553 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
6e4947dd 554 snprintf(key,55,"hMCPhotBoth_cen%d",cent) ;
555 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
556 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
557 snprintf(key,55,"hMCPhotBothcore_cen%d",cent) ;
558 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
559 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
7063ee4d 560 snprintf(key,55,"hMCPhotBoth2_cen%d",cent) ;
561 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
562 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
563 snprintf(key,55,"hMCPhotBoth2core_cen%d",cent) ;
564 fOutputContainer->Add(new TH1F(key,"dN/dpt" ,nPt,ptMin,ptMax));
565 ((TH1F*)fOutputContainer->Last())->Sumw2() ;
566
567 char phiTitle[15]={"TPC"};
568 TH2F * tmp = 0;
569 tmp=new TH2F(Form("hPhotPhi%sAll_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
570 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
571 tmp=new TH2F(Form("hPhotPhi%sAllcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
572 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
573 tmp=new TH2F(Form("hPhotPhi%sDisp_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
574 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
575 tmp=new TH2F(Form("hPhotPhi%sDispcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
576 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
577 tmp=new TH2F(Form("hPhotPhi%sDisp2_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
578 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
579 tmp=new TH2F(Form("hPhotPhi%sDisp2core_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
580 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
581 tmp=new TH2F(Form("hPhotPhi%sCPV_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
582 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
583 tmp=new TH2F(Form("hPhotPhi%sCPVcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
584 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
585 tmp=new TH2F(Form("hPhotPhi%sCPV2_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
586 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
587 tmp=new TH2F(Form("hPhotPhi%sCPV2core_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
588 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
589 tmp=new TH2F(Form("hPhotPhi%sBoth_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
590 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
591 tmp=new TH2F(Form("hPhotPhi%sBothcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
592 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
593 tmp=new TH2F(Form("hPhotPhi%sBoth2_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
594 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
595 tmp=new TH2F(Form("hPhotPhi%sBoth2core_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPtf,xPt,nPhi,xPhi) ;
596 tmp->Sumw2() ; fOutputContainer->Add(tmp) ;
597
02f0b8e7 598
7063ee4d 599 //Pions for flow - with weight 1/Nclu
600 TH3F * tmp3 = 0;
601 tmp3 = new TH3F(Form("hMassPt%sAll_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
602 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
603 tmp3 = new TH3F(Form("hMassPt%sAllcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
604 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
605 tmp3 = new TH3F(Form("hMassPt%sCPV_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
606 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
607 tmp3 = new TH3F(Form("hMassPt%sCPVcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
608 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
609 tmp3 = new TH3F(Form("hMassPt%sCPV2_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
610 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
611 tmp3 = new TH3F(Form("hMassPt%sCPV2core_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
612 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
613 tmp3 = new TH3F(Form("hMassPt%sDisp_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
614 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
615 tmp3 = new TH3F(Form("hMassPt%sDispcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
616 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
617 tmp3 = new TH3F(Form("hMassPt%sDisp2_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
618 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
619 tmp3 = new TH3F(Form("hMassPt%sDisp2core_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
620 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
621 tmp3 = new TH3F(Form("hMassPt%sBoth_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
622 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
623 tmp3 = new TH3F(Form("hMassPt%sBothcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
624 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
625 tmp3 = new TH3F(Form("hMassPt%sBoth2_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
626 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
627 tmp3 = new TH3F(Form("hMassPt%sBoth2core_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
628 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
629
630//mixed
631 tmp3 = new TH3F(Form("hMiMassPt%sAll_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
632 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
633 tmp3 = new TH3F(Form("hMiMassPt%sAllcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
634 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
635 tmp3 = new TH3F(Form("hMiMassPt%sCPV_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
636 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
637 tmp3 = new TH3F(Form("hMiMassPt%sCPVcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
638 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
639 tmp3 = new TH3F(Form("hMiMassPt%sCPV2_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
640 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
641 tmp3 = new TH3F(Form("hMiMassPt%sCPV2core_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
642 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
643 tmp3 = new TH3F(Form("hMiMassPt%sDisp_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
644 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
645 tmp3 = new TH3F(Form("hMiMassPt%sDispcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
646 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
647 tmp3 = new TH3F(Form("hMiMassPt%sDisp2_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
648 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
649 tmp3 = new TH3F(Form("hMiMassPt%sDisp2core_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
650 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
651 tmp3 = new TH3F(Form("hMiMassPt%sBoth_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
652 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
653 tmp3 = new TH3F(Form("hMiMassPt%sBothcore_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
654 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
655 tmp3 = new TH3F(Form("hMiMassPt%sBoth2_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
656 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
657 tmp3 = new TH3F(Form("hMiMassPt%sBoth2core_cen%d",phiTitle,cent),"(M,p_{T},d#phi)_{#gamma#gamma}",nMm,xM,nPtf,xPt,nPhi,xPhi);
658 tmp3->Sumw2() ; fOutputContainer->Add(tmp3) ;
659
02f0b8e7 660 }
661
6e4947dd 662 fOutputContainer->Add(new TH2F("hMCPi0M11","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
663 fOutputContainer->Add(new TH2F("hMCPi0M22","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
664 fOutputContainer->Add(new TH2F("hMCPi0M33","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
665 fOutputContainer->Add(new TH2F("hMCPi0M12","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
666 fOutputContainer->Add(new TH2F("hMCPi0M13","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
667 fOutputContainer->Add(new TH2F("hMCPi0M23","(M,p_{T},d#phi)_{#gamma#gamma}" ,nM,mMin,mMax,nPt,ptMin,ptMax));
02f0b8e7 668
669 //MC
670 for(Int_t cent=0; cent<6; cent++){
671 snprintf(key,55,"hMC_rap_gamma_cen%d",cent) ;
672 fOutputContainer->Add(new TH1F(key,"Rapidity pi0",200,-1.,1.)) ;
673 snprintf(key,55,"hMC_rap_pi0_cen%d",cent) ;
674 fOutputContainer->Add(new TH1F(key,"Rapidity pi0",200,-1.,1.)) ;
675 snprintf(key,55,"hMC_rap_eta_cen%d",cent) ;
676 fOutputContainer->Add(new TH1F("hMC_rap_eta","Rapidity eta",200,-1.,1.)) ;
677 snprintf(key,55,"hMC_phi_gamma_cen%d",cent) ;
678 fOutputContainer->Add(new TH1F(key,"Phi pi0",200,0.,TMath::TwoPi())) ;
679 snprintf(key,55,"hMC_phi_pi0_cen%d",cent) ;
680 fOutputContainer->Add(new TH1F(key,"Phi pi0",200,0.,TMath::TwoPi())) ;
681 snprintf(key,55,"hMC_phi_eta_cen%d",cent) ;
682 fOutputContainer->Add(new TH1F(key,"Phi eta",200,0.,TMath::TwoPi())) ;
683 snprintf(key,55,"hMC_all_gamma_cen%d",cent) ;
684 fOutputContainer->Add(new TH1F(key,"Rapidity photon",250,0.,25.)) ;
685 snprintf(key,55,"hMC_all_pi0_cen%d",cent) ;
686 fOutputContainer->Add(new TH1F(key,"Rapidity pi0",250,0.,25.)) ;
687 snprintf(key,55,"hMC_all_eta_cen%d",cent) ;
688 fOutputContainer->Add(new TH1F(key,"Rapidity eta",250,0.,25.)) ;
689 snprintf(key,55,"hMC_unitEta_gamma_cen%d",cent) ;
690 fOutputContainer->Add(new TH1F(key,"Pt photon",250,0.,25.)) ;
691 snprintf(key,55,"hMC_unitEta_pi0_cen%d",cent) ;
692 fOutputContainer->Add(new TH1F(key,"Rapidity eta",250,0.,25.)) ;
693 snprintf(key,55,"hMC_unitEta_eta_cen%d",cent) ;
694 fOutputContainer->Add(new TH1F(key,"Rapidity eta",250,0.,25.)) ;
695 }
696
697 PostData(1, fOutputContainer);
698
699}
700
701//________________________________________________________________________
702void AliAnalysisTaskPi0DiffEfficiency::UserExec(Option_t *)
703{
704 // Main loop, called for each event
705 // Analyze ESD/AOD
706
707 FillHistogram("hSelEvents",0.5) ;
708
709 AliAODEvent *event = dynamic_cast<AliAODEvent*>(InputEvent());
710 if (!event) {
711 Printf("ERROR: Could not retrieve event");
712 PostData(1, fOutputContainer);
713 return;
714 }
715
716 FillHistogram("hSelEvents",1.5) ;
0a918d8d 717 AliAODHeader *header = dynamic_cast<AliAODHeader*>(event->GetHeader()) ;
718 if(!header) AliFatal("Not a standard AOD");
02f0b8e7 719
720 // Checks if we have a primary vertex
721 // Get primary vertices form ESD
722 const AliAODVertex *esdVertex5 = event->GetPrimaryVertex();
723
724 // don't rely on ESD vertex, assume (0,0,0)
725 Double_t vtx0[3] ={0.,0.,0.};
02f0b8e7 726
727
728 FillHistogram("hZvertex",esdVertex5->GetZ());
729 if (TMath::Abs(esdVertex5->GetZ()) > 10. ){
730 PostData(1, fOutputContainer);
731 return;
732 }
733 FillHistogram("hSelEvents",2.5) ;
734
735 //Vtx class z-bin
736 // Int_t zvtx = (Int_t)((vtx5[2]+10.)/2.) ;
737 // if(zvtx<0)zvtx=0 ;
738 // if(zvtx>9)zvtx=9 ;
739 Int_t zvtx=0 ;
740
741// fCentrality=header->GetCentralityP()->GetCentralityPercentile("V0M"); // returns the centrality percentile,
742// //a float from 0 to 100 (or to the trigger efficiency)
743 fCentrality=header->GetZDCN2Energy() ;
744
745 if( fCentrality < 0. || fCentrality>80.){
746 PostData(1, fOutputContainer);
747 return;
748 }
749 FillHistogram("hSelEvents",3.5) ;
750 Float_t bins[7]={0.,5.,10.,20.,40.,60.,80.} ;
751 fCenBin=0 ;
752 while(fCenBin<6 && fCentrality > bins[fCenBin+1])
753 fCenBin++ ;
754
7063ee4d 755 const Int_t nMixEvents[6]={4,4,5,10,20,20} ;
756
02f0b8e7 757
7063ee4d 758 //reaction plane
02f0b8e7 759 fRPfull= header->GetZDCN1Energy() ;
760 if(fRPfull==999){ //reaction plain was not defined
761 PostData(1, fOutputContainer);
762 return;
763 }
764
765 FillHistogram("hSelEvents",4.5) ;
766 //All event selections done
767 FillHistogram("hCentrality",fCentrality) ;
768 //Reaction plain is defined in the range (-pi/2;pi/2)
769 //We have 10 bins
6e4947dd 770 Int_t irp=Int_t(10.*fRPfull/TMath::Pi());
02f0b8e7 771 if(irp>9)irp=9 ;
772
773 if(!fPHOSEvents[zvtx][fCenBin][irp])
774 fPHOSEvents[zvtx][fCenBin][irp]=new TList() ;
775 TList * prevPHOS = fPHOSEvents[zvtx][fCenBin][irp] ;
776
777 // Get PHOS rotation matrices from ESD and set them to the PHOS geometry
778 if(fEventCounter == 0) {
779 for(Int_t mod=0; mod<5; mod++) {
780 const TGeoHMatrix* m =header->GetPHOSMatrix(mod) ;
781 fPHOSGeo->SetMisalMatrix(m,mod) ;
782 Printf("PHOS geo matrix for module # %d is set: %p\n", mod,m);
783 }
784 fEventCounter++ ;
785 }
6631c482 786
787 TClonesArray *mcArray = (TClonesArray*)event->FindListObject(AliAODMCParticle::StdBranchName());
02f0b8e7 788
789 ProcessMC() ;
790
791 if(fPHOSEvent1){
792 fPHOSEvent1->Clear() ;
793 fPHOSEvent2->Clear() ;
794 }
795 else{
796 fPHOSEvent1 = new TClonesArray("AliCaloPhoton",200) ;
797 fPHOSEvent2 = new TClonesArray("AliCaloPhoton",200) ;
798 }
799
800 TClonesArray * clustersEmb = (TClonesArray*)event->FindListObject("EmbeddedCaloClusters") ;
6e4947dd 801 AliAODCaloCells * cellsEmb = (AliAODCaloCells *)event->FindListObject("EmbeddedPHOScells") ;
02f0b8e7 802 TClonesArray * clustersOld = event->GetCaloClusters() ;
6e4947dd 803 AliAODCaloCells * cellsOld = event->GetPHOSCells() ;
6631c482 804// TClonesArray *mcArray = (TClonesArray*)event->FindListObject(AliAODMCParticle::StdBranchName());
7063ee4d 805
806
02f0b8e7 807 TVector3 vertex(vtx0);
808 char key[55] ;
809 //Before Embedding
810 Int_t multClustOld = clustersOld->GetEntriesFast();
811 Int_t multClustEmb = clustersEmb->GetEntriesFast();
812 Int_t inPHOSold=0 ;
813 Int_t inPHOSemb=0 ;
814 for (Int_t i=0; i<multClustOld; i++) {
815 AliAODCaloCluster *clu = (AliAODCaloCluster*)clustersOld->At(i);
816 if ( !clu->IsPHOS() || clu->E()<0.3) continue;
817
818 Bool_t survive=kFALSE ;
819 for(Int_t ii=0;(ii<multClustEmb)&&(!survive);ii++){
820 AliAODCaloCluster *clu2 = (AliAODCaloCluster*)clustersEmb->At(ii);
821 survive=IsSameCluster(clu,clu2);
822 }
823
824
825 Float_t position[3];
826 clu->GetPosition(position);
827 TVector3 global(position) ;
828 Int_t relId[4] ;
829 fPHOSGeo->GlobalPos2RelId(global,relId) ;
830 Int_t mod = relId[0] ;
831 Int_t cellX = relId[2];
832 Int_t cellZ = relId[3] ;
833 if ( !IsGoodChannel("PHOS",mod,cellX,cellZ) )
834 continue ;
835 if(clu->GetNCells()<3)
836 continue ;
837
1f77ed3b 838 snprintf(key,55,"hCluM%d",mod) ;
02f0b8e7 839 FillHistogram(key,cellX,cellZ,1.);
840
841 TLorentzVector pv1 ;
842 clu->GetMomentum(pv1 ,vtx0);
843
7063ee4d 844 pv1*=Scale(pv1.E()) ;
845
02f0b8e7 846 if(inPHOSold>=fPHOSEvent1->GetSize()){
847 fPHOSEvent1->Expand(inPHOSold+50) ;
848 }
849 AliCaloPhoton * ph = new((*fPHOSEvent1)[inPHOSold]) AliCaloPhoton(pv1.X(),pv1.Py(),pv1.Z(),pv1.E()) ;
850 ph->SetModule(mod) ;
6e4947dd 851 AliPHOSAodCluster cluPHOS1(*clu);
852 cluPHOS1.Recalibrate(fPHOSCalibData,cellsOld); // modify the cell energies
853 Double_t ecore=CoreEnergy(&cluPHOS1) ;
7063ee4d 854 ecore*=Scale(ecore) ;
6e4947dd 855 pv1*= ecore/pv1.E() ;
02f0b8e7 856 ph->SetMomV2(&pv1) ;
857 ph->SetNCells(clu->GetNCells());
7063ee4d 858 Double_t m02=0.,m20=0.;
859 EvalLambdas(&cluPHOS1,0,m02, m20);
860 ph->SetDispBit(TestLambda(clu->E(),m20,m02)) ;
861 EvalLambdas(&cluPHOS1,1,m02, m20);
862 ph->SetDisp2Bit(TestLambda2(clu->E(),m20,m02)) ;
863
6e4947dd 864 ph->SetCPVBit(clu->GetEmcCpvDistance()>2.) ;
865 ph->SetCPV2Bit(clu->GetEmcCpvDistance()>4.) ;
02f0b8e7 866 if(!survive) //this cluster found in list after embedding, skipping it
6e4947dd 867 ph->SetTagged(1) ;
868 ph->SetPhoton(clu->GetNExMax()<2); // Remember, if it is unfolded
7063ee4d 869 ph->SetWeight(1.) ; //All weights for real particles ==1.
02f0b8e7 870
cf150629 871 if(!survive){
872 Double_t distBC=clu->GetDistanceToBadChannel();
873 if(distBC>2.)
874 FillHistogram(Form("hPhotAll_DistBad2_cen%d",fCenBin),ph->Pt(),-1.) ;
875 if(distBC>4.)
876 FillHistogram(Form("hPhotAll_DistBad4_cen%d",fCenBin),ph->Pt(),-1.) ;
877 if(distBC>6.)
878 FillHistogram(Form("hPhotAll_DistBad6_cen%d",fCenBin),ph->Pt(),-1.) ;
879 }
02f0b8e7 880 inPHOSold++ ;
881 }
882
883 for (Int_t i=0; i<multClustEmb; i++) {
884 AliAODCaloCluster *clu = (AliAODCaloCluster*)clustersEmb->At(i);
885 if ( !clu->IsPHOS() || clu->E()<0.3) continue;
886
887 Bool_t survive=kFALSE ;
888 for(Int_t ii=0;(ii<multClustOld)&&(!survive);ii++){
889 AliAODCaloCluster *clu2 = (AliAODCaloCluster*)clustersOld->At(ii);
890 survive=IsSameCluster(clu,clu2);
891 }
892
893 Float_t position[3];
894 clu->GetPosition(position);
895 TVector3 global(position) ;
896 Int_t relId[4] ;
897 fPHOSGeo->GlobalPos2RelId(global,relId) ;
898 Int_t mod = relId[0] ;
899 Int_t cellX = relId[2];
900 Int_t cellZ = relId[3] ;
901 if ( !IsGoodChannel("PHOS",mod,cellX,cellZ) )
902 continue ;
903 if(clu->GetNCells()<3)
904 continue ;
905
1f77ed3b 906 snprintf(key,55,"hCluM%d",mod) ;
02f0b8e7 907 FillHistogram(key,cellX,cellZ,1.);
908
909 TLorentzVector pv1 ;
910 clu->GetMomentum(pv1 ,vtx0);
7063ee4d 911
912 pv1*=Scale(pv1.E()) ;
02f0b8e7 913
914 if(inPHOSemb>=fPHOSEvent2->GetSize()){
915 fPHOSEvent2->Expand(inPHOSemb+50) ;
916 }
917 AliCaloPhoton * ph = new((*fPHOSEvent2)[inPHOSemb]) AliCaloPhoton(pv1.X(),pv1.Py(),pv1.Z(),pv1.E()) ;
918 ph->SetModule(mod) ;
6e4947dd 919 AliPHOSAodCluster cluPHOS1(*clu);
920 cluPHOS1.Recalibrate(fPHOSCalibData,cellsEmb); // modify the cell energies
921 Double_t ecore=CoreEnergy(&cluPHOS1) ;
7063ee4d 922 ecore*=Scale(ecore) ;
6e4947dd 923 pv1*= ecore/pv1.E() ;
02f0b8e7 924 ph->SetMomV2(&pv1) ;
925 ph->SetNCells(clu->GetNCells());
7063ee4d 926 Double_t m02=0.,m20=0.;
927 EvalLambdas(&cluPHOS1,0,m02, m20);
928 ph->SetDispBit(TestLambda(clu->E(),m20,m02)) ;
929 EvalLambdas(&cluPHOS1,1,m02, m20);
930 ph->SetDisp2Bit(TestLambda2(clu->E(),m20,m02)) ;
931
6e4947dd 932 ph->SetCPVBit(clu->GetEmcCpvDistance()>2.) ;
933 ph->SetCPV2Bit(clu->GetEmcCpvDistance()>4.) ;
02f0b8e7 934 if(!survive) //this cluster found in list after embedding, skipping it
6e4947dd 935 ph->SetTagged(1) ;
936 ph->SetPhoton(clu->GetNExMax()<2); // Remember, if it is unfolded
7063ee4d 937
938 //Set weight for embedded particles
939 Double_t w=1. ;
940 Int_t iprim = clu->GetLabel() ;
941 if(iprim<mcArray->GetEntriesFast() && iprim>-1){
942 AliAODMCParticle* particle = (AliAODMCParticle*) mcArray->At(iprim);
943 iprim=particle->GetMother() ;
944 while(iprim>-1){
945 particle = (AliAODMCParticle*) mcArray->At(iprim);
946 iprim=particle->GetMother() ;
947 }
948 if(particle->GetPdgCode()==111){
949 Double_t pt = particle->Pt() ;
950 w=PrimaryWeight(pt) ;
951 }
952 }
953 ph->SetWeight(w) ;
954
02f0b8e7 955
cf150629 956 if(!survive){
957 Double_t distBC=clu->GetDistanceToBadChannel();
958 if(distBC>2.)
959 FillHistogram(Form("hPhotAll_DistBad2_cen%d",fCenBin),ph->Pt()) ;
960 if(distBC>4.)
961 FillHistogram(Form("hPhotAll_DistBad4_cen%d",fCenBin),ph->Pt()) ;
962 if(distBC>6.)
963 FillHistogram(Form("hPhotAll_DistBad6_cen%d",fCenBin),ph->Pt()) ;
964 }
965
02f0b8e7 966 inPHOSemb++ ;
967 }
968
969 //Single photon
970 for (Int_t i1=0; i1<inPHOSold; i1++) {
971 AliCaloPhoton * ph1=(AliCaloPhoton*)fPHOSEvent1->At(i1) ;
972 if(!ph1->IsTagged())
973 continue ;
7063ee4d 974
975 Double_t dphi=ph1->Phi()-fRPfull ;
976 while(dphi<0)dphi+=TMath::Pi() ;
977 while(dphi>TMath::Pi())dphi-=TMath::Pi() ;
978 Double_t pt = ph1->Pt() ;
979 Double_t ptV2=ph1->GetMomV2()->Pt() ;
980//always 1! Double_t w=ph1->GetWeight() ;
981
982 FillHistogram(Form("hPhotPhiTPCAll_cen%d",fCenBin),pt,dphi,-1.) ;
983 FillHistogram(Form("hPhotPhiTPCAllcore_cen%d",fCenBin),ptV2,dphi,-1.) ;
984
985 FillHistogram(Form("hPhotAll_cen%d",fCenBin),pt,-1.) ;
986 FillHistogram(Form("hNegPhotAll_cen%d",fCenBin),pt,-1.) ;
987 FillHistogram(Form("hPhotAllcore_cen%d",fCenBin),ptV2,-1.) ;
988 FillHistogram(Form("hNegPhotAllcore_cen%d",fCenBin),ptV2,-1.) ;
6e4947dd 989 if(ph1->IsPhoton()){
7063ee4d 990 FillHistogram(Form("hPhotAllwou_cen%d",fCenBin),pt,-1.) ;
991 }
02f0b8e7 992 if(ph1->IsCPVOK() ){
7063ee4d 993 FillHistogram(Form("hPhotCPV_cen%d",fCenBin),pt,-1.) ;
994 FillHistogram(Form("hNegPhotCPV_cen%d",fCenBin),pt,-1.) ;
995 FillHistogram(Form("hPhotCPVcore_cen%d",fCenBin),ptV2,-1.) ;
996 FillHistogram(Form("hNegPhotCPVcore_cen%d",fCenBin),ptV2,-1.) ;
997 FillHistogram(Form("hPhotPhiTPCCPV_cen%d",fCenBin),pt,dphi,-1.) ;
998 FillHistogram(Form("hPhotPhiTPCCPVcore_cen%d",fCenBin),ptV2,dphi,-1.) ;
6e4947dd 999 }
1000 if(ph1->IsCPV2OK() ){
7063ee4d 1001 FillHistogram(Form("hPhotCPV2_cen%d",fCenBin),pt,-1.) ;
1002 FillHistogram(Form("hPhotCPV2core_cen%d",fCenBin),ptV2,-1.) ;
1003 FillHistogram(Form("hNegPhotCPV2_cen%d",fCenBin),pt,-1.) ;
1004 FillHistogram(Form("hPhotPhiTPCCPV2_cen%d",fCenBin),pt,dphi,-1.) ;
1005 FillHistogram(Form("hPhotPhiTPCCPV2core_cen%d",fCenBin),ptV2,dphi,-1.) ;
6e4947dd 1006 }
1007 if(ph1->IsDisp2OK()){
7063ee4d 1008 FillHistogram(Form("hPhotDisp2_cen%d",fCenBin),pt,-1.) ;
1009 FillHistogram(Form("hPhotDisp2core_cen%d",fCenBin),ptV2,-1.) ;
1010 FillHistogram(Form("hNegPhotDisp2_cen%d",fCenBin),pt,-1.) ;
1011 FillHistogram(Form("hNegPhotDisp2core_cen%d",fCenBin),ptV2,-1.) ;
1012
1013 FillHistogram(Form("hPhotPhiTPCDisp2_cen%d",fCenBin),pt,dphi,-1.) ;
1014 FillHistogram(Form("hPhotPhiTPCDisp2core_cen%d",fCenBin),ptV2,dphi,-1.) ;
1015 if(ph1->IsCPVOK()){
1016 FillHistogram(Form("hPhotBoth2_cen%d",fCenBin),pt,-1.) ;
1017 FillHistogram(Form("hNegPhotBoth2_cen%d",fCenBin),pt,-1.) ;
1018 FillHistogram(Form("hPhotBoth2core_cen%d",fCenBin),ptV2,-1.) ;
1019 FillHistogram(Form("hNegPhotBoth2core_cen%d",fCenBin),ptV2,-1.) ;
1020 FillHistogram(Form("hPhotPhiTPCBoth2_cen%d",fCenBin),pt,dphi,-1.) ;
1021 FillHistogram(Form("hPhotPhiTPCBoth2core_cen%d",fCenBin),ptV2,dphi,-1.) ;
1022 }
02f0b8e7 1023 }
1024 if(ph1->IsDispOK()){
7063ee4d 1025 FillHistogram(Form("hPhotDisp_cen%d",fCenBin),pt,-1.) ;
1026 FillHistogram(Form("hNegPhotDisp_cen%d",fCenBin),pt,-1.) ;
1027 FillHistogram(Form("hPhotDispcore_cen%d",fCenBin),ptV2,-1.) ;
1028 FillHistogram(Form("hNegPhotDispcore_cen%d",fCenBin),ptV2,-1.) ;
1029 if(ph1->IsPhoton()){
1030 FillHistogram(Form("hPhotDispwou_cen%d",fCenBin),pt,-1.) ;
1031 }
1032
1033 FillHistogram(Form("hPhotPhiTPCDisp_cen%d",fCenBin),pt,dphi,-1.) ;
1034 FillHistogram(Form("hPhotPhiTPCDispcore_cen%d",fCenBin),ptV2,dphi,-1.) ;
02f0b8e7 1035 if(ph1->IsCPVOK()){
7063ee4d 1036 FillHistogram(Form("hPhotBoth_cen%d",fCenBin),pt,-1.) ;
1037 FillHistogram(Form("hNegPhotBoth_cen%d",fCenBin),pt,-1.) ;
1038 FillHistogram(Form("hPhotBothcore_cen%d",fCenBin),ptV2,-1.) ;
1039 FillHistogram(Form("hNegPhotBothcore_cen%d",fCenBin),ptV2,-1.) ;
1040
1041 FillHistogram(Form("hPhotPhiTPCBoth_cen%d",fCenBin),pt,dphi,-1.) ;
1042 FillHistogram(Form("hPhotPhiTPCBothcore_cen%d",fCenBin),ptV2,dphi,-1.) ;
02f0b8e7 1043 }
7063ee4d 1044 }
02f0b8e7 1045 } // end of loop i1
1046
1047 for (Int_t i1=0; i1<inPHOSemb; i1++) {
1048 AliCaloPhoton * ph1=(AliCaloPhoton*)fPHOSEvent2->At(i1) ;
1049 if(!ph1->IsTagged())
1050 continue ;
7063ee4d 1051
1052 Double_t dphi=ph1->Phi()-fRPfull ;
1053 while(dphi<0)dphi+=TMath::Pi() ;
1054 while(dphi>TMath::Pi())dphi-=TMath::Pi() ;
1055 Double_t pt = ph1->Pt() ;
1056 Double_t ptV2=ph1->GetMomV2()->Pt() ;
1057 Double_t w=ph1->GetWeight() ;
1058
1059 FillHistogram(Form("hPhotPhiTPCAll_cen%d",fCenBin),pt,dphi,w) ;
1060 FillHistogram(Form("hPhotPhiTPCAllcore_cen%d",fCenBin),ptV2,dphi,w) ;
1061
1062 FillHistogram(Form("hPhotAll_cen%d",fCenBin),pt,w) ;
1063 FillHistogram(Form("hPhotAllcore_cen%d",fCenBin),ptV2,w) ;
6e4947dd 1064 if(ph1->IsPhoton()){
7063ee4d 1065 FillHistogram(Form("hPhotAllwou_cen%d",fCenBin),pt,w) ;
6e4947dd 1066 }
02f0b8e7 1067 if(ph1->IsCPVOK() ){
7063ee4d 1068 FillHistogram(Form("hPhotCPV_cen%d",fCenBin),pt,w) ;
1069 FillHistogram(Form("hPhotCPVcore_cen%d",fCenBin),ptV2,w) ;
1070
1071 FillHistogram(Form("hPhotPhiTPCCPV_cen%d",fCenBin),pt,dphi,w) ;
1072 FillHistogram(Form("hPhotPhiTPCCPVcore_cen%d",fCenBin),ptV2,dphi,w) ;
6e4947dd 1073 }
1074 if(ph1->IsCPV2OK() ){
7063ee4d 1075 FillHistogram(Form("hPhotCPV2_cen%d",fCenBin),pt,w) ;
1076 FillHistogram(Form("hPhotCPV2core_cen%d",fCenBin),ptV2,w) ;
1077 FillHistogram(Form("hPhotPhiTPCCPV2_cen%d",fCenBin),pt,dphi,w) ;
1078 FillHistogram(Form("hPhotPhiTPCCPV2core_cen%d",fCenBin),ptV2,dphi,w) ;
6e4947dd 1079 }
1080 if(ph1->IsDisp2OK()){
7063ee4d 1081 FillHistogram(Form("hPhotDisp2_cen%d",fCenBin),pt,w) ;
1082 FillHistogram(Form("hPhotDisp2core_cen%d",fCenBin),ptV2,w) ;
1083
1084 FillHistogram(Form("hPhotPhiTPCDisp2_cen%d",fCenBin),pt,dphi,w) ;
1085 FillHistogram(Form("hPhotPhiTPCDisp2core_cen%d",fCenBin),ptV2,dphi,w) ;
1086 if(ph1->IsCPVOK()){
1087 FillHistogram(Form("hPhotBoth2_cen%d",fCenBin),pt,w) ;
1088 FillHistogram(Form("hPhotBoth2core_cen%d",fCenBin),ptV2,w) ;
1089
1090 FillHistogram(Form("hPhotPhiTPCBoth2_cen%d",fCenBin),pt,dphi,w) ;
1091 FillHistogram(Form("hPhotPhiTPCBoth2core_cen%d",fCenBin),ptV2,dphi,w) ;
1092 }
02f0b8e7 1093 }
1094 if(ph1->IsDispOK()){
7063ee4d 1095 FillHistogram(Form("hPhotDisp_cen%d",fCenBin),pt,w) ;
1096 FillHistogram(Form("hPhotDispcore_cen%d",fCenBin),ptV2,w) ;
1097 if(ph1->IsPhoton()){
1098 FillHistogram(Form("hPhotDispwou_cen%d",fCenBin),pt,w) ;
1099 }
1100
1101 FillHistogram(Form("hPhotPhiTPCBoth_cen%d",fCenBin),pt,dphi,w) ;
1102 FillHistogram(Form("hPhotPhiTPCBothcore_cen%d",fCenBin),ptV2,dphi,w) ;
02f0b8e7 1103 if(ph1->IsCPVOK()){
7063ee4d 1104 FillHistogram(Form("hPhotBoth_cen%d",fCenBin),pt,w) ;
1105 FillHistogram(Form("hPhotBothcore_cen%d",fCenBin),ptV2,w) ;
1106
1107 FillHistogram(Form("hPhotPhiTPCBoth_cen%d",fCenBin),pt,dphi,w) ;
1108 FillHistogram(Form("hPhotPhiTPCBothcore_cen%d",fCenBin),ptV2,dphi,w) ;
02f0b8e7 1109 }
1110 } // end of loop i2
1111 } // end of loop i1
1112
1113
7063ee4d 1114 const Double_t prob[10]={0.1,0.2,0.3,1.,1.,1.,1.,1.,1.,1.} ; //Probabilities to accept Tagged+Bg pair
1115
02f0b8e7 1116
1117 // Fill Real disribution:
1118 // Disappeared clusters enter with negative contribution
1119 // In addition fill control histogam with Real before embedding
1120 for (Int_t i1=0; i1<inPHOSold-1; i1++) {
1121 AliCaloPhoton * ph1=(AliCaloPhoton*)fPHOSEvent1->At(i1) ;
7063ee4d 1122 Double_t w1 = ph1->GetWeight() ;
02f0b8e7 1123 for (Int_t i2=i1+1; i2<inPHOSold; i2++) {
1124 AliCaloPhoton * ph2=(AliCaloPhoton*)fPHOSEvent1->At(i2) ;
1125
1126 TLorentzVector p12 = *ph1 + *ph2;
1127 TLorentzVector pv12 = *(ph1->GetMomV2()) + *(ph2->GetMomV2());
7063ee4d 1128 Double_t w2 = ph2->GetWeight() ;
1129 Double_t w = TMath::Sqrt(w1*w2) ;
1130
1131 Double_t dphi=p12.Phi()-fRPfull ;
1132 while(dphi<0)dphi+=TMath::Pi() ;
1133 while(dphi>TMath::Pi())dphi-=TMath::Pi() ;
1134 Double_t m=p12.M() ;
1135 Double_t mV2=pv12.M() ;
1136 Double_t pt = p12.Pt() ;
1137 Double_t ptV2 = pv12.Pt() ;
1138
1139
02f0b8e7 1140 //Fill Controll histogram: Real before embedding
7063ee4d 1141 FillHistogram(Form("hOldMassPtAll_cen%d",fCenBin),m,pt,-w) ;
1142 FillHistogram(Form("hOldMassPtAllcore_cen%d",fCenBin),mV2,ptV2,-w) ;
02f0b8e7 1143 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
7063ee4d 1144 FillHistogram(Form("hOldMassPtCPV_cen%d",fCenBin),m ,pt,-w) ;
1145 FillHistogram(Form("hOldMassPtCPV_cen%d",fCenBin),mV2, ptV2,-w) ;
6e4947dd 1146 }
1147 if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
7063ee4d 1148 FillHistogram(Form("hOldMassPtCPV2_cen%d",fCenBin),m ,pt,-w) ;
6e4947dd 1149 }
1150 if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
7063ee4d 1151 FillHistogram(Form("hOldMassPtDisp2_cen%d",fCenBin),m ,pt,-w) ;
1152 FillHistogram(Form("hOldMassPtDisp2core_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1153 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1154 FillHistogram(Form("hOldMassPtBoth2_cen%d",fCenBin),m ,pt,-w) ;
1155 FillHistogram(Form("hOldMassPtBoth2core_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1156 }
02f0b8e7 1157 }
1158 if(ph1->IsDispOK() && ph2->IsDispOK()){
7063ee4d 1159 FillHistogram(Form("hOldMassPtDisp_cen%d",fCenBin),m ,pt,-w) ;
1160 FillHistogram(Form("hOldMassPtDispcore_cen%d",fCenBin),mV2 ,ptV2,-w) ;
02f0b8e7 1161 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
7063ee4d 1162 FillHistogram(Form("hOldMassPtBoth_cen%d",fCenBin),m ,pt,-w) ;
1163 FillHistogram(Form("hOldMassPtBothcore_cen%d",fCenBin),mV2 ,ptV2,-w) ;
02f0b8e7 1164 }
1165 }
1166
6e4947dd 1167 //Now fill main histograms with negative contributions
02f0b8e7 1168 if(!(ph1->IsTagged() || ph2->IsTagged()) )
1169 continue ;
7063ee4d 1170 if(!ph1->IsTagged() || !ph2->IsTagged()){ //Tagged + Bg combination
1171 if(gRandom->Uniform()>prob[fCenBin])
1172 continue ;
6e4947dd 1173 }
7063ee4d 1174 FillHistogram(Form("hMassPtAll_cen%d",fCenBin),m ,pt,-w) ;
1175 FillHistogram(Form("hMassPtAllcore_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1176 if(ph1->IsPhoton()&&ph2->IsPhoton()){
1177 FillHistogram(Form("hMassPtAllwou_cen%d",fCenBin),m,pt,-w) ;
6e4947dd 1178 }
7063ee4d 1179
1180 FillHistogram(Form("hMassPtTPCAll_cen%d",fCenBin),m,pt,dphi,-w) ;
1181 FillHistogram(Form("hMassPtTPCAllcore_cen%d",fCenBin),mV2,ptV2,dphi,-w) ;
1182
1183 FillHistogram(Form("hNegMassPtAll_cen%d",fCenBin),m ,pt,-w) ;
1184 FillHistogram(Form("hNegMassPtAllcore_cen%d",fCenBin),mV2 ,ptV2,-w) ;
6e4947dd 1185
02f0b8e7 1186 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
7063ee4d 1187 FillHistogram(Form("hMassPtCPV_cen%d",fCenBin),m ,pt,-w) ;
1188 FillHistogram(Form("hMassPtCPVcore_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1189
1190 FillHistogram(Form("hMassPtTPCCPV_cen%d",fCenBin),m,pt,dphi,-w) ;
1191 FillHistogram(Form("hMassPtTPCCPVcore_cen%d",fCenBin),mV2,ptV2,dphi,-w) ;
1192
1193 FillHistogram(Form("hNegMassPtCPV_cen%d",fCenBin),m ,pt,-w) ;
1194 FillHistogram(Form("hNegMassPtCPVcore_cen%d",fCenBin),mV2 ,ptV2,-w) ;
6e4947dd 1195 }
1196 if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
7063ee4d 1197 FillHistogram(Form("hMassPtCPV2_cen%d",fCenBin),m ,pt,-w) ;
1198 FillHistogram(Form("hNegMassPtCPV2_cen%d",fCenBin),m ,pt,-w) ;
1199
1200 FillHistogram(Form("hMassPtTPCCPV2_cen%d",fCenBin),m,pt,dphi,-w) ;
1201 FillHistogram(Form("hMassPtTPCCPV2core_cen%d",fCenBin),mV2,ptV2,dphi,-w) ;
6e4947dd 1202 }
1203
1204 if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
7063ee4d 1205 FillHistogram(Form("hMassPtDisp2_cen%d",fCenBin),m ,pt,-w) ;
1206 FillHistogram(Form("hMassPtDisp2core_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1207
1208 FillHistogram(Form("hMassPtTPCDisp2_cen%d",fCenBin),m,pt,dphi,-w) ;
1209 FillHistogram(Form("hMassPtTPCDisp2core_cen%d",fCenBin),mV2,ptV2,dphi,-w) ;
1210
1211 FillHistogram(Form("hNegMassPtDisp2_cen%d",fCenBin),m ,pt,-w) ;
1212 FillHistogram(Form("hNegMassPtDisp2core_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1213 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1214 FillHistogram(Form("hMassPtBoth2_cen%d",fCenBin),m ,pt,-w) ;
1215 FillHistogram(Form("hMassPtBoth2core_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1216 FillHistogram(Form("hNegMassPtBoth2_cen%d",fCenBin),m ,pt,-w) ;
1217 FillHistogram(Form("hNegMassPtBoth2core_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1218
1219 FillHistogram(Form("hMassPtTPCBoth2_cen%d",fCenBin),m,pt,dphi,-w) ;
1220 FillHistogram(Form("hMassPtTPCBoth2core_cen%d",fCenBin),mV2,ptV2,dphi,-w) ;
1221 }
02f0b8e7 1222 }
1223 if(ph1->IsDispOK() && ph2->IsDispOK()){
7063ee4d 1224 FillHistogram(Form("hMassPtDisp_cen%d",fCenBin),m ,pt,-w) ;
1225 FillHistogram(Form("hMassPtDispcore_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1226 if(ph1->IsPhoton()&&ph2->IsPhoton()){
1227 FillHistogram(Form("hMassPtDispwou_cen%d",fCenBin),m,pt,-w) ;
6e4947dd 1228 }
7063ee4d 1229 FillHistogram(Form("hNegMassPtDisp_cen%d",fCenBin),m ,pt,-w) ;
1230 FillHistogram(Form("hNegMassPtDispcore_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1231
1232 FillHistogram(Form("hMassPtTPCDisp_cen%d",fCenBin),m,pt,dphi,-w) ;
1233 FillHistogram(Form("hMassPtTPCDispcore_cen%d",fCenBin),mV2,ptV2,dphi,-w) ;
02f0b8e7 1234 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
7063ee4d 1235 FillHistogram(Form("hMassPtBoth_cen%d",fCenBin),m ,pt,-w) ;
1236 FillHistogram(Form("hMassPtBothcore_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1237 FillHistogram(Form("hNegMassPtBoth_cen%d",fCenBin),m ,pt,-w) ;
1238 FillHistogram(Form("hNegMassPtBothcore_cen%d",fCenBin),mV2 ,ptV2,-w) ;
1239
1240 FillHistogram(Form("hMassPtTPCBoth_cen%d",fCenBin),m,pt,dphi,-w) ;
1241 FillHistogram(Form("hMassPtTPCBothcore_cen%d",fCenBin),mV2,ptV2,dphi,-w) ;
02f0b8e7 1242 }
1243 }
1244 } // end of loop i2
1245 } // end of loop i1
1246
1247
1248 // Further fill Real disribution
1249 // now with positive contribution from new clusters
1250 // ass well fill controll histogram
1251 for (Int_t i1=0; i1<inPHOSemb-1; i1++) {
1252 AliCaloPhoton * ph1=(AliCaloPhoton*)fPHOSEvent2->At(i1) ;
7063ee4d 1253 Double_t w1 = ph1->GetWeight() ;
02f0b8e7 1254 for (Int_t i2=i1+1; i2<inPHOSemb; i2++) {
1255 AliCaloPhoton * ph2=(AliCaloPhoton*)fPHOSEvent2->At(i2) ;
1256
1257 TLorentzVector p12 = *ph1 + *ph2;
1258 TLorentzVector pv12 = *(ph1->GetMomV2()) + *(ph2->GetMomV2());
7063ee4d 1259 Double_t m=p12.M() ;
1260 Double_t mV2=pv12.M() ;
1261 Double_t pt = p12.Pt() ;
1262 Double_t ptV2 = pv12.Pt() ;
1263 Double_t w2 = ph2->GetWeight() ;
1264 Double_t w = TMath::Sqrt(w1*w2) ;
02f0b8e7 1265
1266 // Controll histogram: Real after embedding
7063ee4d 1267 FillHistogram(Form("hNewMassPtAll_cen%d",fCenBin),m ,pt,w) ;
1268 FillHistogram(Form("hNewMassPtAllcore_cen%d",fCenBin),mV2 ,ptV2,w) ;
02f0b8e7 1269 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
7063ee4d 1270 FillHistogram(Form("hNewMassPtCPV_cen%d",fCenBin),m ,pt,w) ;
1271 FillHistogram(Form("hNewMassPtCPVcore_cen%d",fCenBin),mV2 ,ptV2,w) ;
6e4947dd 1272 }
1273 if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
7063ee4d 1274 FillHistogram(Form("hNewMassPtCPV2_cen%d",fCenBin),m ,pt,w) ;
6e4947dd 1275 }
1276 if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
7063ee4d 1277 FillHistogram(Form("hNewMassPtDisp2_cen%d",fCenBin),m ,pt,w) ;
1278 FillHistogram(Form("hNewMassPtDisp2core_cen%d",fCenBin),mV2 ,ptV2,w) ;
1279 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1280 FillHistogram(Form("hNewMassPtBoth2_cen%d",fCenBin),m ,pt,w) ;
1281 FillHistogram(Form("hNewMassPtBoth2core_cen%d",fCenBin),mV2,ptV2,w) ;
1282 }
02f0b8e7 1283 }
1284 if(ph1->IsDispOK() && ph2->IsDispOK()){
7063ee4d 1285 FillHistogram(Form("hNewMassPtDisp_cen%d",fCenBin),m ,pt,w) ;
1286 FillHistogram(Form("hNewMassPtDispcore_cen%d",fCenBin),mV2 ,ptV2,w) ;
02f0b8e7 1287 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
7063ee4d 1288 FillHistogram(Form("hNewMassPtBoth_cen%d",fCenBin),m ,pt,w) ;
1289 FillHistogram(Form("hNewMassPtBothcore_cen%d",fCenBin),mV2 ,ptV2,w) ;
02f0b8e7 1290 }
1291 }
1292
1293 //Now fill main histogamm
1294 //new clusters with positive contribution
1295 if(!(ph1->IsTagged() || ph2->IsTagged()) )
1296 continue ;
7063ee4d 1297 if(!ph1->IsTagged() || !ph2->IsTagged()){ //Tagged + Bg combination
1298 if(gRandom->Uniform()>prob[fCenBin])
1299 continue ;
6e4947dd 1300 }
7063ee4d 1301
1302 Double_t dphi=p12.Phi()-fRPfull ;
1303 while(dphi<0)dphi+=TMath::Pi() ;
1304 while(dphi>TMath::Pi())dphi-=TMath::Pi() ;
1305
1306 FillHistogram(Form("hMassPtAll_cen%d",fCenBin),m ,pt,w) ;
1307 FillHistogram(Form("hMassPtAllcore_cen%d",fCenBin),mV2 ,ptV2,w) ;
1308 if(ph1->IsPhoton()&&ph2->IsPhoton()){
1309 FillHistogram(Form("hMassPtAllwou_cen%d",fCenBin),m,pt,w) ;
6e4947dd 1310 }
7063ee4d 1311
1312 FillHistogram(Form("hMassPtTPCAll_cen%d",fCenBin),m,pt,dphi,w) ;
1313 FillHistogram(Form("hMassPtTPCAllcore_cen%d",fCenBin),mV2,ptV2,dphi,w) ;
1314
02f0b8e7 1315 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
7063ee4d 1316 FillHistogram(Form("hMassPtCPV_cen%d",fCenBin),m ,pt,w) ;
1317 FillHistogram(Form("hMassPtCPVcore_cen%d",fCenBin),mV2 ,ptV2,w) ;
1318
1319 FillHistogram(Form("hMassPtTPCCPV_cen%d",fCenBin),m,pt,dphi,w) ;
1320 FillHistogram(Form("hMassPtTPCCPVcore_cen%d",fCenBin),mV2,ptV2,dphi,w) ;
6e4947dd 1321 }
1322 if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
7063ee4d 1323 FillHistogram(Form("hMassPtCPV2_cen%d",fCenBin),m ,pt,w) ;
1324
1325 FillHistogram(Form("hMassPtTPCCPV2_cen%d",fCenBin),m,pt,dphi,w) ;
1326 FillHistogram(Form("hMassPtTPCCPV2core_cen%d",fCenBin),mV2,ptV2,dphi,w) ;
6e4947dd 1327 }
1328 if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
7063ee4d 1329 FillHistogram(Form("hMassPtDisp2_cen%d",fCenBin),m ,pt,w) ;
1330 FillHistogram(Form("hMassPtDisp2core_cen%d",fCenBin),mV2 ,ptV2,w) ;
1331
1332 FillHistogram(Form("hMassPtTPCDisp2_cen%d",fCenBin),m,pt,dphi,w) ;
1333 FillHistogram(Form("hMassPtTPCDisp2core_cen%d",fCenBin),mV2,ptV2,dphi,w) ;
1334 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1335 FillHistogram(Form("hMassPtBoth2_cen%d",fCenBin),m ,pt,w) ;
1336 FillHistogram(Form("hMassPtBoth2core_cen%d",fCenBin),mV2 ,ptV2,w) ;
1337
1338 FillHistogram(Form("hMassPtTPCBoth2_cen%d",fCenBin),m,pt,dphi,w) ;
1339 FillHistogram(Form("hMassPtTPCBoth2core_cen%d",fCenBin),mV2,ptV2,dphi,w) ;
1340 }
1341 }
02f0b8e7 1342 if(ph1->IsDispOK() && ph2->IsDispOK()){
7063ee4d 1343 FillHistogram(Form("hMassPtDisp_cen%d",fCenBin),m ,pt,w) ;
1344 FillHistogram(Form("hMassPtDispcore_cen%d",fCenBin),mV2 ,ptV2,w) ;
1345 if(ph1->IsPhoton()&&ph2->IsPhoton()){
1346 FillHistogram(Form("hMassPtDispwou_cen%d",fCenBin),m,pt,w) ;
6e4947dd 1347 }
02f0b8e7 1348
7063ee4d 1349 FillHistogram(Form("hMassPtTPCDisp_cen%d",fCenBin),m,pt,dphi,w) ;
1350 FillHistogram(Form("hMassPtTPCDispcore_cen%d",fCenBin),mV2,ptV2,dphi,w) ;
1351
02f0b8e7 1352 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
7063ee4d 1353 FillHistogram(Form("hMassPtBoth_cen%d",fCenBin),m ,pt,w) ;
1354 FillHistogram(Form("hMassPtBothcore_cen%d",fCenBin),mV2 ,ptV2,w) ;
1355
1356 FillHistogram(Form("hMassPtTPCBoth_cen%d",fCenBin),m,pt,dphi,w) ;
1357 FillHistogram(Form("hMassPtTPCBothcore_cen%d",fCenBin),mV2,ptV2,dphi,w) ;
02f0b8e7 1358 }
1359 }
1360 } // end of loop i2
1361 } // end of loop i1
1362
1363
1364 //now mixed, does not really matter old or new list of clusters
1365 for (Int_t i1=0; i1<inPHOSemb; i1++) {
1366 AliCaloPhoton * ph1=(AliCaloPhoton*)fPHOSEvent2->At(i1) ;
1367 for(Int_t ev=0; ev<prevPHOS->GetSize();ev++){
1368 TClonesArray * mixPHOS = static_cast<TClonesArray*>(prevPHOS->At(ev)) ;
1369 for(Int_t i2=0; i2<mixPHOS->GetEntriesFast();i2++){
1370 AliCaloPhoton * ph2=(AliCaloPhoton*)mixPHOS->At(i2) ;
1371 TLorentzVector p12 = *ph1 + *ph2;
1372 TLorentzVector pv12 = *(ph1->GetMomV2()) + *(ph2->GetMomV2());
7063ee4d 1373
1374 Double_t dphi=p12.Phi()-fRPfull ;
1375 while(dphi<0)dphi+=TMath::Pi() ;
1376 while(dphi>TMath::Pi())dphi-=TMath::Pi() ;
1377 Double_t m=p12.M() ;
1378 Double_t mV2=pv12.M() ;
1379 Double_t pt = p12.Pt() ;
1380 Double_t ptV2 = pv12.Pt() ;
02f0b8e7 1381
7063ee4d 1382 FillHistogram(Form("hMiMassPtAll_cen%d",fCenBin),m ,pt,1.) ;
1383 FillHistogram(Form("hMiMassPtAllcore_cen%d",fCenBin),mV2 ,ptV2,1.) ;
1384 if(ph1->IsPhoton()&&ph2->IsPhoton()){
1385 FillHistogram(Form("hMiMassPtAllwou_cen%d",fCenBin),m,pt,1.) ;
6e4947dd 1386 }
7063ee4d 1387
1388 FillHistogram(Form("hMiMassPtTPCAll_cen%d",fCenBin),m,pt,dphi) ;
1389 FillHistogram(Form("hMiMassPtTPCAllcore_cen%d",fCenBin),mV2,ptV2,dphi) ;
1390
02f0b8e7 1391 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
7063ee4d 1392 FillHistogram(Form("hMiMassPtCPV_cen%d",fCenBin),m ,pt,1.) ;
1393 FillHistogram(Form("hMiMassPtCPVcore_cen%d",fCenBin),mV2 ,ptV2,1.) ;
1394
1395 FillHistogram(Form("hMiMassPtTPCCPV_cen%d",fCenBin),m,pt,dphi) ;
1396 FillHistogram(Form("hMiMassPtTPCCPVcore_cen%d",fCenBin),mV2,ptV2,dphi) ;
1397
6e4947dd 1398 }
1399 if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
7063ee4d 1400 FillHistogram(Form("hMiMassPtCPV2_cen%d",fCenBin),m ,pt,1.) ;
1401
1402 FillHistogram(Form("hMiMassPtTPCCPV2_cen%d",fCenBin),m,pt,dphi) ;
1403 FillHistogram(Form("hMiMassPtTPCCPV2core_cen%d",fCenBin),mV2,ptV2,dphi) ;
6e4947dd 1404 }
1405 if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
7063ee4d 1406 FillHistogram(Form("hMiMassPtDisp2_cen%d",fCenBin),m ,pt,1.) ;
1407 FillHistogram(Form("hMiMassPtDisp2core_cen%d",fCenBin),mV2 ,ptV2,1.) ;
1408
1409 FillHistogram(Form("hMiMassPtTPCDisp2_cen%d",fCenBin),m,pt,dphi) ;
1410 FillHistogram(Form("hMiMassPtTPCDisp2core_cen%d",fCenBin),mV2,ptV2,dphi) ;
1411
1412 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1413 FillHistogram(Form("hMiMassPtBoth2_cen%d",fCenBin),m ,pt,1.) ;
1414 FillHistogram(Form("hMiMassPtBoth2core_cen%d",fCenBin),mV2 ,ptV2,1.) ;
1415
1416 FillHistogram(Form("hMiMassPtTPCBoth2_cen%d",fCenBin),m,pt,dphi) ;
1417 FillHistogram(Form("hMiMassPtTPCBoth2core_cen%d",fCenBin),mV2,ptV2,dphi) ;
1418 }
02f0b8e7 1419 }
1420 if(ph1->IsDispOK() && ph2->IsDispOK()){
7063ee4d 1421 FillHistogram(Form("hMiMassPtDisp_cen%d",fCenBin),m ,pt,1.) ;
1422 FillHistogram(Form("hMiMassPtDispcore_cen%d",fCenBin),mV2 ,ptV2,1.) ;
1423 if(ph1->IsPhoton()&&ph2->IsPhoton()){
1424 FillHistogram(Form("hMiMassPtDispwou_cen%d",fCenBin),m,pt,1.) ;
1425 }
1426
1427 FillHistogram(Form("hMiMassPtTPCDisp_cen%d",fCenBin),m,pt,dphi) ;
1428 FillHistogram(Form("hMiMassPtTPCDispcore_cen%d",fCenBin),mV2,ptV2,dphi) ;
02f0b8e7 1429
1430 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
7063ee4d 1431 FillHistogram(Form("hMiMassPtBoth_cen%d",fCenBin),m ,pt,1.) ;
1432 FillHistogram(Form("hMiMassPtBothcore_cen%d",fCenBin),mV2 ,ptV2,1.) ;
1433
1434 FillHistogram(Form("hMiMassPtTPCBoth_cen%d",fCenBin),m,pt,dphi) ;
1435 FillHistogram(Form("hMiMassPtTPCBothcore_cen%d",fCenBin),mV2,ptV2,dphi) ;
02f0b8e7 1436 }
1437 }
1438 } // end of loop i2
1439 }
1440 } // end of loop i1
7063ee4d 1441
02f0b8e7 1442
1443
1444 //Now we either add current events to stack or remove
1445 //If no photons in current event - no need to add it to mixed
1446 if(fPHOSEvent2->GetEntriesFast()>0){
1447 prevPHOS->AddFirst(fPHOSEvent2) ;
1448 fPHOSEvent2=0;
1449 delete fPHOSEvent1;
1450 fPHOSEvent1=0;
7063ee4d 1451 if(prevPHOS->GetSize()>nMixEvents[fCenBin]){//Remove redundant events
02f0b8e7 1452 TClonesArray * tmp = static_cast<TClonesArray*>(prevPHOS->Last()) ;
1453 prevPHOS->RemoveLast() ;
1454 delete tmp ;
1455 }
1456 }
1457 // Post output data.
1458 PostData(1, fOutputContainer);
1459 fEventCounter++;
1460}
02f0b8e7 1461//___________________________________________________________________________
1462Bool_t AliAnalysisTaskPi0DiffEfficiency::IsSameCluster(AliAODCaloCluster * c1, AliAODCaloCluster * c2)const{
1463 //Compare clusters before and after embedding
1464 //clusters are the same if
1465 // - Energy changed less than 0.1% (numerical accuracy in reconstruction)
1466 // - lists of digits are the same
1467
1468 if(c1->GetNCells() != c2->GetNCells())
1469 return kFALSE ;
1470
6e4947dd 1471 if(TMath::Abs(c1->E()-c2->E())>0.01*c1->E())
02f0b8e7 1472 return kFALSE ;
1473
1474 UShort_t *list1 = c1->GetCellsAbsId() ;
1475 UShort_t *list2 = c2->GetCellsAbsId() ;
1476 for(Int_t i=0; i<c1->GetNCells(); i++){
1477 if(list1[i] != list2[i])
1478 return kFALSE ;
1479 }
1480 return kTRUE ;
1481
1482}
7063ee4d 1483//____________________________________________________________________________
1484void AliAnalysisTaskPi0DiffEfficiency::EvalLambdas(AliAODCaloCluster * clu, Int_t iR,Double_t &m02, Double_t &m20){
1485 //calculate dispecrsion of the cluster in the circle with radius distanceCut around the maximum
1486
1487 Double_t rCut=0. ;
1488 if(iR==0)
1489 rCut=3.5 ;
1490 else
1491 rCut=4.5 ;
1492
1493
1494 Double32_t * elist = clu->GetCellsAmplitudeFraction() ;
1495// Calculates the center of gravity in the local PHOS-module coordinates
1496 Float_t wtot = 0;
1497 Double_t xc[100]={0} ;
1498 Double_t zc[100]={0} ;
1499 Double_t x = 0 ;
1500 Double_t z = 0 ;
1501 Int_t mulDigit=TMath::Min(100,clu->GetNCells()) ;
1502 const Double_t logWeight=4.5 ;
1503 for(Int_t iDigit=0; iDigit<mulDigit; iDigit++) {
1504 Int_t relid[4] ;
1505 Float_t xi ;
1506 Float_t zi ;
1507 fPHOSGeo->AbsToRelNumbering(clu->GetCellAbsId(iDigit), relid) ;
1508 fPHOSGeo->RelPosInModule(relid, xi, zi);
1509 xc[iDigit]=xi ;
1510 zc[iDigit]=zi ;
1511 if (clu->E()>0 && elist[iDigit]>0) {
1512 Float_t w = TMath::Max( 0., logWeight + TMath::Log( elist[iDigit] / clu->E() ) ) ;
1513 x += xc[iDigit] * w ;
1514 z += zc[iDigit] * w ;
1515 wtot += w ;
1516 }
1517 }
1518 if (wtot>0) {
1519 x /= wtot ;
1520 z /= wtot ;
1521 }
1522
1523 wtot = 0. ;
1524 Double_t dxx = 0.;
1525 Double_t dzz = 0.;
1526 Double_t dxz = 0.;
1527 Double_t xCut = 0. ;
1528 Double_t zCut = 0. ;
1529 for(Int_t iDigit=0; iDigit<mulDigit; iDigit++) {
1530 if (clu->E()>0 && elist[iDigit]>0.) {
1531 Double_t w = TMath::Max( 0., logWeight + TMath::Log( elist[iDigit] / clu->E() ) ) ;
1532 Double_t xi= xc[iDigit] ;
1533 Double_t zi= zc[iDigit] ;
1534 if((xi-x)*(xi-x)+(zi-z)*(zi-z) < rCut*rCut){
1535 xCut += w * xi ;
1536 zCut += w * zi ;
1537 dxx += w * xi * xi ;
1538 dzz += w * zi * zi ;
1539 dxz += w * xi * zi ;
1540 wtot += w ;
1541 }
1542 }
1543
1544 }
1545 if (wtot>0) {
1546 xCut/= wtot ;
1547 zCut/= wtot ;
1548 dxx /= wtot ;
1549 dzz /= wtot ;
1550 dxz /= wtot ;
1551 dxx -= xCut * xCut ;
1552 dzz -= zCut * zCut ;
1553 dxz -= xCut * zCut ;
1554
1555 m02 = 0.5 * (dxx + dzz) + TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz ) ;
1556 m20 = 0.5 * (dxx + dzz) - TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz ) ;
1557 }
1558 else {
1559 m20=m02=0.;
1560 }
1561
1562}
1563//___________________________________________________________________________
1564void AliAnalysisTaskPi0DiffEfficiency::ProcessMC(){
1565 //fill histograms for efficiensy etc. calculation
1566 const Double_t rcut = 1. ; //cut for primary particles
1567 //---------First pi0/eta-----------------------------
1568 char partName[10] ;
1569
1570 AliAODEvent *event = dynamic_cast<AliAODEvent*>(InputEvent());
1571 if(!event) return ;
1572 TClonesArray *mcArray = (TClonesArray*)event->FindListObject(AliAODMCParticle::StdBranchName());
1573 for(Int_t i=0;i<mcArray->GetEntriesFast();i++){
1574 AliAODMCParticle* particle = (AliAODMCParticle*) mcArray->At(i);
1575 if(particle->GetPdgCode() == 111)
1576 snprintf(partName,10,"pi0") ;
1577 else
1578 if(particle->GetPdgCode() == 221)
1579 snprintf(partName,10,"eta") ;
1580 else
1581 if(particle->GetPdgCode() == 22)
1582 snprintf(partName,10,"gamma") ;
1583 else
1584 continue ;
1585
1586 //Primary particle
1587 Double_t r=TMath::Sqrt(particle->Xv()*particle->Xv()+particle->Yv()*particle->Yv());
1588 if(r >rcut)
1589 continue ;
1590
1591 Double_t pt = particle->Pt() ;
1592 Double_t w = PrimaryWeight(pt) ;
1593 //Total number of pi0 with creation radius <1 cm
1594 FillHistogram(Form("hMC_all_%s_cen%d",partName,fCenBin),pt,w) ;
1595 if(TMath::Abs(particle->Y())<0.12){
1596 FillHistogram(Form("hMC_unitEta_%s_cen%d",partName,fCenBin),pt,w) ;
1597 }
1598
1599 FillHistogram(Form("hMC_rap_%s_cen%d",partName,fCenBin),particle->Y(),w) ;
1600
1601 Double_t phi=particle->Phi() ;
1602 while(phi<0.)phi+=TMath::TwoPi() ;
1603 while(phi>TMath::TwoPi())phi-=TMath::TwoPi() ;
1604 FillHistogram(Form("hMC_phi_%s_cen%d",partName,fCenBin),phi,w) ;
1605
1606 }
1607
1608 //Now calculate "Real" distribution of clusters with primary
1609 TClonesArray * clusters = (TClonesArray*)event->FindListObject("EmbeddedCaloClusters") ;
1610 AliAODCaloCells * cellsEmb = (AliAODCaloCells *)event->FindListObject("EmbeddedPHOScells") ;
1611 Int_t multClust = clusters->GetEntriesFast();
1612 TClonesArray cluPrim("AliCaloPhoton",multClust) ; //clusters with primary
1613 Int_t inPHOS=0 ;
1614 Double_t vtx0[3] = {0,0,0};
1615 for (Int_t i=0; i<multClust; i++) {
1616 AliAODCaloCluster *clu = (AliAODCaloCluster*)clusters->At(i);
1617 if ( !clu->IsPHOS() || clu->E()<0.3) continue;
1618 if(clu->GetLabel()<0) continue ;
1619
1620 Float_t position[3];
1621 clu->GetPosition(position);
1622 TVector3 global(position) ;
1623 Int_t relId[4] ;
1624 fPHOSGeo->GlobalPos2RelId(global,relId) ;
1625 Int_t mod = relId[0] ;
1626 Int_t cellX = relId[2];
1627 Int_t cellZ = relId[3] ;
1628 if ( !IsGoodChannel("PHOS",mod,cellX,cellZ) )
1629 continue ;
1630 if(clu->GetNCells()<3)
1631 continue ;
1632
1633 TLorentzVector pv1 ;
1634 clu->GetMomentum(pv1 ,vtx0);
1635
1636 pv1*=Scale(pv1.E()) ;
1637
1638
1639
1640 if(inPHOS>=cluPrim.GetSize()){
1641 cluPrim.Expand(inPHOS+50) ;
1642 }
1643 AliCaloPhoton * ph = new(cluPrim[inPHOS]) AliCaloPhoton(pv1.X(),pv1.Py(),pv1.Z(),pv1.E()) ;
1644 //AliCaloPhoton * ph = (AliCaloPhoton*)fPHOSEvent->At(inPHOS) ;
1645 ph->SetModule(mod) ;
1646
1647 AliPHOSAodCluster cluPHOS1(*clu);
1648 cluPHOS1.Recalibrate(fPHOSCalibData,cellsEmb); // modify the cell energies
1649 Double_t ecore=CoreEnergy(&cluPHOS1) ;
1650 ecore*=Scale(ecore) ;
1651 pv1*= ecore/pv1.E() ;
1652 ph->SetMomV2(&pv1) ;
1653 ph->SetNCells(clu->GetNCells());
1654 Double_t m02=0.,m20=0.;
1655 EvalLambdas(&cluPHOS1,0,m02, m20);
1656 ph->SetDispBit(TestLambda(clu->E(),m20,m02)) ;
1657 EvalLambdas(&cluPHOS1,1,m02, m20);
1658 ph->SetDisp2Bit(TestLambda2(clu->E(),m20,m02)) ;
1659
1660 ph->SetCPVBit(clu->GetEmcCpvDistance()>2.) ; //radius in sigmas
1661 ph->SetCPV2Bit(clu->GetEmcCpvDistance()>4.) ;
1662 ph->SetPhoton(clu->GetNExMax()<2); // Remember, if it is unfolded
1663 Double_t w=1. ;
1664 Int_t iprim = clu->GetLabel() ;
1665 if(iprim<mcArray->GetEntriesFast() && iprim>-1){
1666 AliAODMCParticle* particle = (AliAODMCParticle*) mcArray->At(iprim);
1667 iprim=particle->GetMother() ;
1668 while(iprim>-1){
1669 particle = (AliAODMCParticle*) mcArray->At(iprim);
1670 iprim=particle->GetMother() ;
1671 }
1672 if(particle->GetPdgCode()==111){
1673 Double_t pt = particle->Pt() ;
1674 w=PrimaryWeight(pt) ;
1675 }
1676 }
1677 ph->SetWeight(w) ;
1678
1679 inPHOS++ ;
1680 }
1681
1682 //Single photon
1683 char key[55] ;
1684 for (Int_t i1=0; i1<inPHOS; i1++) {
1685 AliCaloPhoton * ph1=(AliCaloPhoton*)cluPrim.At(i1) ;
1686 Double_t pt=ph1->Pt() ;
1687 Double_t ptV2=ph1->GetMomV2()->Pt() ;
1688 Double_t w = ph1->GetWeight() ;
1689 FillHistogram(Form("hMCPhotAll_cen%d",fCenBin),pt,w) ;
1690 FillHistogram(Form("hMCPhotAllcore_cen%d",fCenBin),ptV2,w) ;
1691 if(ph1->IsPhoton()){
1692 FillHistogram(Form("hMCPhotAllwou_cen%d",fCenBin),pt,w) ;
1693 }
1694 if(ph1->IsCPVOK() ){
1695 FillHistogram(Form("hMCPhotCPV_cen%d",fCenBin),pt,w) ;
1696 FillHistogram(Form("hMCPhotCPVcore_cen%d",fCenBin),ptV2,w) ;
1697 }
1698 if(ph1->IsCPV2OK() ){
1699 snprintf(key,55,"hMCPhotCPV2_cen%d",fCenBin) ;
1700 FillHistogram(Form("hMCPhotCPV2_cen%d",fCenBin),pt,w) ;
1701 }
1702 if(ph1->IsDisp2OK()){
1703 FillHistogram(Form("hMCPhotDisp2_cen%d",fCenBin),pt,w) ;
1704 FillHistogram(Form("hMCPhotDisp2core_cen%d",fCenBin),ptV2,w) ;
1705 if(ph1->IsCPVOK()){
1706 FillHistogram(Form("hMCPhotBoth2_cen%d",fCenBin),pt,w) ;
1707 FillHistogram(Form("hMCPhotBoth2core_cen%d",fCenBin),ptV2,w) ;
1708 }
1709 }
1710 if(ph1->IsDispOK()){
1711 FillHistogram(Form("hMCPhotDisp_cen%d",fCenBin),pt,w) ;
1712 FillHistogram(Form("hMCPhotDispcore_cen%d",fCenBin),ptV2,w) ;
1713 if(ph1->IsPhoton()){
1714 FillHistogram(Form("hMCPhotDispwou_cen%d",fCenBin),pt,w) ;
1715 }
1716 if(ph1->IsCPVOK()){
1717 FillHistogram(Form("hMCPhotBoth_cen%d",fCenBin),pt,w) ;
1718 FillHistogram(Form("hMCPhotBothcore_cen%d",fCenBin),ptV2,w) ;
1719 }
1720 } // end of loop i2
1721 } // end of loop i1
1722
1723 // Fill Real disribution
1724 for (Int_t i1=0; i1<inPHOS-1; i1++) {
1725 AliCaloPhoton * ph1=(AliCaloPhoton*)cluPrim.At(i1) ;
1726 Double_t w1 = ph1->GetWeight() ;
1727 for (Int_t i2=i1+1; i2<inPHOS; i2++) {
1728 AliCaloPhoton * ph2=(AliCaloPhoton*)cluPrim.At(i2) ;
1729 Double_t w2 = ph2->GetWeight() ;
1730 Double_t w = TMath::Sqrt(w1*w2) ;
1731 TLorentzVector p12 = *ph1 + *ph2;
1732 TLorentzVector pv12 = *(ph1->GetMomV2()) + *(ph2->GetMomV2());
1733 Double_t m=p12.M() ;
1734 Double_t pt=p12.Pt() ;
1735 Double_t mV2=pv12.M() ;
1736 Double_t ptV2=pv12.Pt() ;
1737
1738 FillHistogram(Form("hMCMassPtAll_cen%d",fCenBin),m,pt,w) ;
1739 snprintf(key,55,"hMCMassPtAllcore_cen%d",fCenBin) ;
1740 FillHistogram(Form("hMCMassPtAllcore_cen%d",fCenBin),mV2,ptV2,w) ;
1741 if(ph1->IsPhoton()&&ph2->IsPhoton() ){
1742 FillHistogram(Form("hMCMassPtAllwou_cen%d",fCenBin),m,pt,w) ;
1743 }
1744
1745 if(ph1->Module()==1 && ph2->Module()==1)
1746 FillHistogram("hMCPi0M11",m,pt,w);
1747 else if(ph1->Module()==2 && ph2->Module()==2)
1748 FillHistogram("hMCPi0M22",m,pt,w);
1749 else if(ph1->Module()==3 && ph2->Module()==3)
1750 FillHistogram("hMCPi0M33",m,pt,w);
1751 else if(ph1->Module()==1 && ph2->Module()==2)
1752 FillHistogram("hMCPi0M12",m,pt,w);
1753 else if(ph1->Module()==1 && ph2->Module()==3)
1754 FillHistogram("hMCPi0M13",m,pt,w);
1755 else if(ph1->Module()==2 && ph2->Module()==3)
1756 FillHistogram("hMCPi0M23",m,pt,w);
1757
1758 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1759 FillHistogram(Form("hMCMassPtCPV_cen%d",fCenBin),m,pt,w) ;
1760 FillHistogram(Form("hMCMassPtCPVcore_cen%d",fCenBin),mV2,ptV2,w) ;
1761 }
1762 if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
1763 FillHistogram(Form("hMCMassPtCPV2core_cen%d",fCenBin),mV2,ptV2,w) ;
1764 }
1765 if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
1766 FillHistogram(Form("hMCMassPtDisp2_cen%d",fCenBin),m,pt,w) ;
1767 FillHistogram(Form("hMCMassPtDisp2core_cen%d",fCenBin),mV2,ptV2,w) ;
1768 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1769 FillHistogram(Form("hMCMassPtBoth2_cen%d",fCenBin),m,pt,w) ;
1770 FillHistogram(Form("hMCMassPtBoth2core_cen%d",fCenBin),mV2,ptV2,w) ;
1771 }
1772 }
1773 if(ph1->IsDispOK() && ph2->IsDispOK()){
1774 FillHistogram(Form("hMCMassPtDisp_cen%d",fCenBin),m,pt,w) ;
1775 FillHistogram(Form("hMCMassPtDispcore_cen%d",fCenBin),mV2,ptV2,w) ;
1776 if(ph1->IsPhoton()&& ph2->IsPhoton()){
1777 FillHistogram(Form("hMCMassPtDispwou_cen%d",fCenBin),m,pt,w) ;
1778 }
1779 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
1780 FillHistogram(Form("hMCMassPtBoth_cen%d",fCenBin),m,pt,w) ;
1781 FillHistogram(Form("hMCMassPtBothcore_cen%d",fCenBin),mV2,ptV2,w) ;
1782 }
1783 }
1784 } // end of loop i2
1785 } // end of loop i1
1786}
1787//____________________________________________________________________________
1788Double_t AliAnalysisTaskPi0DiffEfficiency::CoreEnergy(AliPHOSAodCluster * clu){
1789 //calculate energy of the cluster in the circle with radius distanceCut around the maximum
1790
1791 //Can not use already calculated coordinates?
1792 //They have incidence correction...
1793 const Double_t distanceCut =3.5 ;
1794 const Double_t logWeight=4.5 ;
1795
1796 Double32_t * elist = clu->GetCellsAmplitudeFraction() ;
1797// Calculates the center of gravity in the local PHOS-module coordinates
1798 Float_t wtot = 0;
1799 Double_t xc[100]={0} ;
1800 Double_t zc[100]={0} ;
1801 Double_t x = 0 ;
1802 Double_t z = 0 ;
1803 Int_t mulDigit=TMath::Min(100,clu->GetNCells()) ;
1804 for(Int_t iDigit=0; iDigit<mulDigit; iDigit++) {
1805 Int_t relid[4] ;
1806 Float_t xi ;
1807 Float_t zi ;
1808 fPHOSGeo->AbsToRelNumbering(clu->GetCellAbsId(iDigit), relid) ;
1809 fPHOSGeo->RelPosInModule(relid, xi, zi);
1810 xc[iDigit]=xi ;
1811 zc[iDigit]=zi ;
1812 if (clu->E()>0 && elist[iDigit]>0) {
1813 Float_t w = TMath::Max( 0., logWeight + TMath::Log( elist[iDigit] / clu->E() ) ) ;
1814 x += xc[iDigit] * w ;
1815 z += zc[iDigit] * w ;
1816 wtot += w ;
1817 }
1818 }
1819 if (wtot>0) {
1820 x /= wtot ;
1821 z /= wtot ;
1822 }
1823 Double_t coreE=0. ;
1824 for(Int_t iDigit=0; iDigit < mulDigit; iDigit++) {
1825 Double_t distance = TMath::Sqrt((xc[iDigit]-x)*(xc[iDigit]-x)+(zc[iDigit]-z)*(zc[iDigit]-z)) ;
1826 if(distance < distanceCut)
1827 coreE += elist[iDigit] ;
1828 }
1829 //Apply non-linearity correction
1830 return (0.0241+1.0504*coreE+0.000249*coreE*coreE) ;
1831}
1832//_____________________________________________________________________________
1833Bool_t AliAnalysisTaskPi0DiffEfficiency::TestLambda(Double_t pt,Double_t l1,Double_t l2){
1834
1835 //For R=3.5
1836 Double_t l1Mean = 1.170014 -0.059465/(1.+0.019343*pt+0.147544*pt*pt) ;
1837 Double_t l2Mean = 1.626270 + 0.761554*exp(-1.213839*pt)-0.020027*pt ;
1838 Double_t l1Sigma = 0.133409 + 0.261307*exp(-0.636874*pt)-0.002849*pt ;
1839 Double_t l2Sigma = 0.289698 + 0.459400*exp(-1.214242*pt)-0.012578*pt ;
1840 Double_t c=-0.124103 ;
1841/*
1842 Double_t l2Mean = 1.53126+9.50835e+06/(1.+1.08728e+07*pt+1.73420e+06*pt*pt) ;
1843 Double_t l1Mean = 1.12365+0.123770*TMath::Exp(-pt*0.246551)+5.30000e-03*pt ;
1844 Double_t l2Sigma = 6.48260e-02+7.60261e+10/(1.+1.53012e+11*pt+5.01265e+05*pt*pt)+9.00000e-03*pt;
1845 Double_t l1Sigma = 4.44719e-04+6.99839e-01/(1.+1.22497e+00*pt+6.78604e-07*pt*pt)+9.00000e-03*pt;
1846 Double_t c=-0.35-0.550*TMath::Exp(-0.390730*pt) ;
1847*/
1848 Double_t R2=0.5*(l1-l1Mean)*(l1-l1Mean)/l1Sigma/l1Sigma +
1849 0.5*(l2-l2Mean)*(l2-l2Mean)/l2Sigma/l2Sigma +
1850 0.5*c*(l1-l1Mean)*(l2-l2Mean)/l1Sigma/l2Sigma ;
1851 return (R2<2.5*2.5) ;
1852
1853}
1854//_____________________________________________________________________________
1855Bool_t AliAnalysisTaskPi0DiffEfficiency::TestLambda2(Double_t pt,Double_t l1,Double_t l2){
1856
1857//For R=4.5
1858 Double_t l1Mean = 1.150200 + 0.097886/(1.+1.486645*pt+0.000038*pt*pt) ;
1859 Double_t l2Mean = 1.574706 + 0.997966*exp(-0.895075*pt)-0.010666*pt ;
1860 Double_t l1Sigma = 0.100255 + 0.337177*exp(-0.517684*pt)+0.001170*pt ;
1861 Double_t l2Sigma = 0.232580 + 0.573401*exp(-0.735903*pt)-0.002325*pt ;
1862 Double_t c = -0.110983 -0.017353/(1.-1.836995*pt+0.934517*pt*pt) ;
1863
1864/*
1865 Double_t l2Mean = 1.53126+9.50835e+06/(1.+1.08728e+07*pt+1.73420e+06*pt*pt) ;
1866 Double_t l1Mean = 1.12365+0.123770*TMath::Exp(-pt*0.246551)+5.30000e-03*pt ;
1867 Double_t l2Sigma = 6.48260e-02+7.60261e+10/(1.+1.53012e+11*pt+5.01265e+05*pt*pt)+9.00000e-03*pt;
1868 Double_t l1Sigma = 4.44719e-04+6.99839e-01/(1.+1.22497e+00*pt+6.78604e-07*pt*pt)+9.00000e-03*pt;
1869 Double_t c=-0.35-0.550*TMath::Exp(-0.390730*pt) ;
1870*/
1871 Double_t R2=0.5*(l1-l1Mean)*(l1-l1Mean)/l1Sigma/l1Sigma +
1872 0.5*(l2-l2Mean)*(l2-l2Mean)/l2Sigma/l2Sigma +
1873 0.5*c*(l1-l1Mean)*(l2-l2Mean)/l1Sigma/l2Sigma ;
1874 return (R2<2.5*2.5) ;
1875
1876}
1877//____________________________________________________________________________
1878Double_t AliAnalysisTaskPi0DiffEfficiency::TestCPV(Double_t dx, Double_t dz, Double_t pt, Int_t charge){
1879 //Parameterization of LHC10h period
1880 //_true if neutral_
1881
1882 Double_t meanX=0;
1883 Double_t meanZ=0.;
1884 Double_t sx=TMath::Min(5.4,2.59719e+02*TMath::Exp(-pt/1.02053e-01)+
1885 6.58365e-01*5.91917e-01*5.91917e-01/((pt-9.61306e-01)*(pt-9.61306e-01)+5.91917e-01*5.91917e-01)+1.59219);
1886 Double_t sz=TMath::Min(2.75,4.90341e+02*1.91456e-02*1.91456e-02/(pt*pt+1.91456e-02*1.91456e-02)+1.60) ;
1887 AliAODEvent *event = dynamic_cast<AliAODEvent*>(InputEvent());
1207640c 1888 if(!event)AliFatal("Can not get ESD event") ;
7063ee4d 1889 Double_t mf = event->GetMagneticField(); //Positive for ++ and negative for --
1890
1891 if(mf<0.){ //field --
1892 meanZ = -0.468318 ;
1893 if(charge>0)
1894 meanX=TMath::Min(7.3, 3.89994*1.20679*1.20679/(pt*pt+1.20679*1.20679)+0.249029+2.49088e+07*TMath::Exp(-pt*3.33650e+01)) ;
1895 else
1896 meanX=-TMath::Min(7.7,3.86040*0.912499*0.912499/(pt*pt+0.912499*0.912499)+1.23114+4.48277e+05*TMath::Exp(-pt*2.57070e+01)) ;
1897 }
1898 else{ //Field ++
1899 meanZ= -0.468318;
1900 if(charge>0)
1901 meanX=-TMath::Min(8.0,3.86040*1.31357*1.31357/(pt*pt+1.31357*1.31357)+0.880579+7.56199e+06*TMath::Exp(-pt*3.08451e+01)) ;
1902 else
1903 meanX= TMath::Min(6.85, 3.89994*1.16240*1.16240/(pt*pt+1.16240*1.16240)-0.120787+2.20275e+05*TMath::Exp(-pt*2.40913e+01)) ;
1904 }
02f0b8e7 1905
7063ee4d 1906 Double_t rz=(dz-meanZ)/sz ;
1907 Double_t rx=(dx-meanX)/sx ;
1908 return TMath::Sqrt(rx*rx+rz*rz) ;
1909}
1910//____________________________________________________________________________
1911Bool_t AliAnalysisTaskPi0DiffEfficiency::TestTOF(Double_t t, Double_t e){
1912
1913 Double_t sigma = TMath::Sqrt(2.23183e-09*2.23183e-09
1914 +2.24611e-09*2.24611e-09/e
1915 +5.65054e-09*5.65054e-09/e/e) ;
1916 sigma=TMath::Min(20.e-9,sigma) ; //for the soft (<400 MeV) part
1917 Double_t mean=1.1e-9 ;
1918 if(TMath::Abs(t-mean)<2.*sigma)
1919 return kTRUE ;
1920 else
1921 if(TMath::Abs(t-mean+100.e-9)<2.*sigma)
1922 return kTRUE ;
1923
1924 return kFALSE ;
1925}
1926//_____________________________________________________________________________
1927void AliAnalysisTaskPi0DiffEfficiency::FillHistogram(const char * key,Double_t x)const{
1928 //FillHistogram
1929 TObject * tmp = fOutputContainer->FindObject(key) ;
1930 if(!tmp){
1931 AliInfo(Form("can not find histogram <%s> ",key)) ;
1932 return ;
1933 }
1934 if(tmp->IsA() == TClass::GetClass("TH1I")){
1935 ((TH1I*)tmp)->Fill(x) ;
1936 return ;
1937 }
1938 if(tmp->IsA() == TClass::GetClass("TH1F")){
1939 ((TH1F*)tmp)->Fill(x) ;
1940 return ;
1941 }
1942 if(tmp->IsA() == TClass::GetClass("TH1D")){
1943 ((TH1D*)tmp)->Fill(x) ;
1944 return ;
1945 }
1946 AliInfo(Form("can not find 1D histogram <%s> ",key)) ;
1947}
1948//_____________________________________________________________________________
1949void AliAnalysisTaskPi0DiffEfficiency::FillHistogram(const char * key,Double_t x,Double_t y)const{
1950 //FillHistogram
1951 TObject * tmp = fOutputContainer->FindObject(key) ;
1952 if(!tmp){
1953 AliInfo(Form("can not find histogram <%s> ",key)) ;
1954 return ;
1955 }
1956 if(tmp->IsA() == TClass::GetClass("TH1F")){
1957 ((TH1F*)tmp)->Fill(x,y) ;
1958 return ;
1959 }
1960 if(tmp->IsA() == TClass::GetClass("TH2F")){
1961 ((TH2F*)tmp)->Fill(x,y) ;
1962 return ;
1963 }
1964 AliError(Form("Calling FillHistogram with 2 parameters for histo <%s> of type %s",key,tmp->IsA()->GetName())) ;
1965}
02f0b8e7 1966
7063ee4d 1967//_____________________________________________________________________________
1968void AliAnalysisTaskPi0DiffEfficiency::FillHistogram(const char * key,Double_t x,Double_t y, Double_t z) const{
1969 //Fills 1D histograms with key
1970 TObject * tmp = fOutputContainer->FindObject(key) ;
1971 if(!tmp){
1972 AliInfo(Form("can not find histogram <%s> ",key)) ;
1973 return ;
1974 }
1975 if(tmp->IsA() == TClass::GetClass("TH2F")){
1976 ((TH2F*)tmp)->Fill(x,y,z) ;
1977 return ;
1978 }
1979 if(tmp->IsA() == TClass::GetClass("TH3F")){
1980 ((TH3F*)tmp)->Fill(x,y,z) ;
1981 return ;
1982 }
1983}
1984//_____________________________________________________________________________
1985void AliAnalysisTaskPi0DiffEfficiency::FillHistogram(const char * key,Double_t x,Double_t y, Double_t z, Double_t w) const{
1986 //Fills 1D histograms with key
1987 TObject * tmp = fOutputContainer->FindObject(key) ;
1988 if(!tmp){
1989 AliInfo(Form("can not find histogram <%s> ",key)) ;
1990 return ;
1991 }
1992 if(tmp->IsA() == TClass::GetClass("TH3F")){
1993 ((TH3F*)tmp)->Fill(x,y,z,w) ;
1994 return ;
1995 }
1996}
1997//_____________________________________________________________________________
1998Bool_t AliAnalysisTaskPi0DiffEfficiency::IsGoodChannel(const char * det, Int_t mod, Int_t ix, Int_t iz)
1999{
2000 //Check if this channel belogs to the good ones
2001
2002 if(strcmp(det,"PHOS")==0){
2003 if(mod>5 || mod<1){
2004 AliError(Form("No bad map for PHOS module %d ",mod)) ;
2005 return kTRUE ;
2006 }
2007 if(!fPHOSBadMap[mod]){
2008 AliError(Form("No Bad map for PHOS module %d",mod)) ;
2009 return kTRUE ;
2010 }
2011 if(fPHOSBadMap[mod]->GetBinContent(ix,iz)>0)
2012 return kFALSE ;
2013 else
2014 return kTRUE ;
2015 }
2016 else{
2017 AliError(Form("Can not find bad channels for detector %s ",det)) ;
2018 }
2019 return kTRUE ;
2020}
2021//_____________________________________________________________________________
2022Double_t AliAnalysisTaskPi0DiffEfficiency::PrimaryWeight(Double_t x){
2023
2024 Double_t w=1 ;
6631c482 2025
2026
2027 //Parameterization of LHC10h data from Jan 2013 (pi0 spectrum)
2028 //Should be consistend with spectrum parameterization used in simulation
7063ee4d 2029 if(fCenBin==0) //0-5
2030 w = (0.561741+0.332841*x-0.007082*x*x)/(1.-0.447804*x+0.157830*x*x)+0.080394*x ;
2031 if(fCenBin==1) //5-10
2032 w = (0.659096+0.101701*x+0.042395*x*x)/(1.-0.470110*x+0.154665*x*x)+0.052932*x ;
2033 if(fCenBin==2) //10-20
2034 w = (0.615575+0.005621*x+0.069263*x*x)/(1.-0.485422*x+0.160822*x*x)+0.040865*x ;
2035 if(fCenBin==3) //20-40
2036 w = (0.441240+0.158358*x+0.059458*x*x)/(1.-0.332609*x+0.147528*x*x)+0.037926*x ;
2037 if(fCenBin==4) //40-60
2038 w = (0.467895-0.001113*x+0.029610*x*x)/(1.-0.266502*x+0.065105*x*x)+0.025431*x ;
2039 if(fCenBin==5) //60-80
2040 w = (0.465204-0.139362*x+0.043500*x*x)/(1.-0.371689*x+0.067606*x*x)+0.006519*x ;
6631c482 2041
2042/*
2043 //Parameterization of photon spectrum 25.02
2044 if(fCenBin==0) //0-5
2045 w=(0.870487-0.494032*x+0.076334*x*x+0.001065*x*x*x)/(1.-0.646014*x+0.113839*x*x);
2046 if(fCenBin==1) //5-10
2047 w=(-8.310403+15.767226*x-2.167756*x*x+0.184356*x*x*x)/(1.+4.556793*x+0.980941*x*x);
2048 if(fCenBin==2) //10-20
2049 w=(-5.281594+7.477165*x-0.688609*x*x+0.097601*x*x*x)/(1.+1.102693*x+0.882454*x*x);
2050 if(fCenBin==3) //20-40
2051 w=(0.789472-4.750155*x+4.381545*x*x-0.029239*x*x*x)/(1.-3.738304*x+3.328318*x*x);
2052 if(fCenBin==4) //40-60
2053 w=(0.876792-0.491379*x+0.130872*x*x-0.001390*x*x*x)/(1.+-0.511934*x+0.112705*x*x);
2054 if(fCenBin==5) //60-80
2055 w=(0.719912-0.167292*x+0.017196*x*x+0.000861*x*x*x)/(1.-0.336146*x+0.037731*x*x);
2056*/
7063ee4d 2057 return TMath::Max(0.,w) ;
2058}
2059
2060
02f0b8e7 2061