bugs corrected
[u/mrichter/AliRoot.git] / PWG4 / PartCorr / AliAnalysisTaskGammaConversion.cxx
CommitLineData
9fdb2372 1/**************************************************************************\r
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
3 * *\r
4 * Author: Ana Marin, Kathrin Koch, Kenneth Aamodt *\r
5 * Version 1.0 *\r
6 * *\r
7 * Permission to use, copy, modify and distribute this software and its *\r
8 * documentation strictly for non-commercial purposes is hereby granted *\r
9 * without fee, provided that the above copyright notice appears in all *\r
10 * copies and that both the copyright notice and this permission notice *\r
11 * appear in the supporting documentation. The authors make no claims *\r
12 * about the suitability of this software for any purpose. It is *\r
13 * provided "as is" without express or implied warranty. *\r
14 **************************************************************************/\r
15\r
16////////////////////////////////////////////////\r
17//--------------------------------------------- \r
18// Class used to do analysis on conversion pairs\r
19//---------------------------------------------\r
20////////////////////////////////////////////////\r
21\r
22// root\r
23#include <TChain.h>\r
24\r
25// analysis\r
26#include "AliAnalysisTaskGammaConversion.h"\r
27#include "AliStack.h"\r
28#include "AliLog.h"\r
29#include <vector>\r
30\r
31class AliKFVertex;\r
32class AliAODHandler;\r
33class AliAODEvent;\r
34class ALiESDEvent;\r
35class AliMCEvent;\r
36class AliMCEventHandler;\r
37class AliESDInputHandler;\r
38class AliAnalysisManager;\r
39class Riostream;\r
40class TFile;\r
41class TInterpreter;\r
42class TSystem;\r
43class TROOT;\r
44\r
45ClassImp(AliAnalysisTaskGammaConversion)\r
46\r
47\r
48AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion():\r
49 AliAnalysisTaskSE(),\r
50 fV0Reader(NULL),\r
51 fStack(NULL),\r
52 fOutputContainer(NULL),\r
53 fHistograms(NULL),\r
54 fDoMCTruth(kFALSE),\r
55 fMCAllGammas(),\r
56 fMCPi0s(),\r
57 fMCEtas(),\r
58 fMCGammaChic(),\r
59 fKFReconstructedGammas(),\r
60 fElectronMass(-1),\r
61 fGammaMass(-1),\r
62 fPi0Mass(-1),\r
63 fEtaMass(-1),\r
64 fGammaWidth(-1),\r
65 fPi0Width(-1),\r
66 fEtaWidth(-1),\r
67 fCalculateBackground(kFALSE)\r
68{\r
69 // Default constructor\r
70 // Common I/O in slot 0\r
71 DefineInput (0, TChain::Class());\r
72 DefineOutput(0, TTree::Class());\r
73\r
74 // Your private output\r
75 DefineOutput(1, TList::Class());\r
76}\r
77\r
78AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion(const char* name):\r
79 AliAnalysisTaskSE(name),\r
80 fV0Reader(NULL),\r
81 fStack(NULL),\r
82 fOutputContainer(0x0),\r
83 fHistograms(NULL),\r
84 fDoMCTruth(kFALSE),\r
85 fMCAllGammas(),\r
86 fMCPi0s(),\r
87 fMCEtas(),\r
88 fMCGammaChic(),\r
89 fKFReconstructedGammas(),\r
90 fElectronMass(-1),\r
91 fGammaMass(-1),\r
92 fPi0Mass(-1),\r
93 fEtaMass(-1),\r
94 fGammaWidth(-1),\r
95 fPi0Width(-1),\r
96 fEtaWidth(-1),\r
97 fCalculateBackground(kFALSE)\r
98{\r
99 // Common I/O in slot 0\r
100 DefineInput (0, TChain::Class());\r
101 DefineOutput(0, TTree::Class());\r
102 \r
103 // Your private output\r
104 DefineOutput(1, TList::Class());\r
105}\r
106\r
107AliAnalysisTaskGammaConversion::~AliAnalysisTaskGammaConversion() \r
108{\r
109 // Remove all pointers\r
110 \r
111 if(fOutputContainer){\r
112 fOutputContainer->Clear() ; \r
113 delete fOutputContainer ;\r
114 }\r
115 if(fHistograms){\r
116 delete fHistograms;\r
117 }\r
118 if(fV0Reader){\r
119 delete fV0Reader;\r
120 }\r
121}\r
122\r
123\r
124void AliAnalysisTaskGammaConversion::Init()\r
125{\r
126 // Initialization\r
a6f44de5 127 AliLog::SetGlobalLogLevel(AliLog::kError);\r
9fdb2372 128}\r
129\r
130\r
131void AliAnalysisTaskGammaConversion::Exec(Option_t */*option*/)\r
132{\r
133 // Execute analysis for current event\r
134 \r
135 ConnectInputData("");\r
136 \r
137 //clear vectors\r
138 fMCAllGammas.clear();\r
139 fMCPi0s.clear();\r
140 fMCEtas.clear();\r
141 fMCGammaChic.clear();\r
142\r
143 fKFReconstructedGammas.clear();\r
144\r
145 //Clear the data in the v0Reader\r
146 fV0Reader->UpdateEventByEventData();\r
147\r
148 // Process the MC information\r
149 if(fDoMCTruth){\r
150 ProcessMCData();\r
151 }\r
152\r
153 // Process the v0 information\r
154 ProcessV0s();\r
155\r
156 //calculate background if flag is set\r
157 if(fCalculateBackground){\r
158 CalculateBackground();\r
159 }\r
160\r
161 // Process reconstructed gammas\r
162 ProcessGammasForNeutralMesonAnalysis();\r
163\r
164 PostData(1, fOutputContainer);\r
165 \r
166}\r
167\r
168void AliAnalysisTaskGammaConversion::ConnectInputData(Option_t */*option*/){\r
169 // see header file for documentation\r
170\r
171 if(fV0Reader == NULL){\r
172 // Write warning here cuts and so on are default if this ever happens\r
173 }\r
174 fV0Reader->Initialize();\r
175}\r
176\r
177void AliAnalysisTaskGammaConversion::ProcessMCData(){\r
178 // see header file for documentation\r
179 \r
180 fStack = fV0Reader->GetMCStack();\r
181\r
182 for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) {\r
183 TParticle* particle = (TParticle *)fStack->Particle(iTracks);\r
184\r
185 if (!particle) {\r
186 //print warning here\r
187 continue;\r
188 }\r
189 \r
190 if(particle->Pt()<fV0Reader->GetPtCut()){\r
191 continue;\r
192 }\r
193 \r
194 if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut()){\r
195 continue;\r
196 }\r
197\r
198 if(particle->R()>fV0Reader->GetMaxRCut()){ // cuts on distance from collision point\r
199 continue;\r
200 }\r
201\r
202 Double_t tmpPhi=particle->Phi();\r
203 if(particle->Phi()> TMath::Pi()){\r
204 tmpPhi = particle->Phi()-(2*TMath::Pi());\r
205 }\r
206\r
207 \r
208 //process the gammas\r
209 if (particle->GetPdgCode()== 22){\r
210 fMCAllGammas.push_back(particle);\r
211 if(particle->GetMother(0)>-1){ //Means we have a mother\r
212 if( fStack->Particle(particle->GetMother(0))->GetPdgCode() != 22 ){//Checks for a non gamma mother.\r
213 fHistograms->FillHistogram("MC_Gamma_Energy", particle->Energy());\r
214 fHistograms->FillHistogram("MC_Gamma_Pt", particle->Pt());\r
215 fHistograms->FillHistogram("MC_Gamma_Eta", particle->Eta());\r
216 \r
217 fHistograms->FillHistogram("MC_Gamma_Phi", tmpPhi);\r
218\r
219 //adding the conversion points from all gammas with e+e- daughters\r
a6f44de5 220 if(particle->GetNDaughters() >= 2){\r
221 TParticle* daughter0 = NULL;\r
222 TParticle* daughter1 = NULL;\r
9fdb2372 223 \r
a6f44de5 224 for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){\r
225 TParticle *tmpDaughter = fStack->Particle(daughterIndex);\r
226 if(tmpDaughter->GetUniqueID() == 5){\r
227 if(tmpDaughter->GetPdgCode() == 11){\r
228 daughter0 = tmpDaughter;\r
229 }\r
230 else if(tmpDaughter->GetPdgCode() == -11){\r
231 daughter1 = tmpDaughter;\r
232 }\r
233 }\r
234 }\r
235\r
236 if(daughter0 == NULL || daughter1 == NULL){ // means we do not have two daughters from pair production\r
237 continue;\r
238 }\r
239\r
9fdb2372 240 if(daughter0->R()>fV0Reader->GetMaxRCut() || daughter1->R()>fV0Reader->GetMaxRCut()){\r
241 continue;\r
242 }\r
243 \r
244 if((daughter0->GetPdgCode() == -11 && daughter1->GetPdgCode()) == 11 ||\r
245 (daughter0->GetPdgCode() == 11 && daughter1->GetPdgCode() == -11)){\r
246\r
247 // begin Mapping \r
248 Int_t rBin = fHistograms->GetRBin(daughter0->R());\r
249 Int_t phiBin = fHistograms->GetPhiBin(daughter0->Phi());\r
250 \r
251 TString nameMCMappingPhiR="";\r
252 nameMCMappingPhiR.Form("MC_EP_Mapping-Phi%02d-R%02d",phiBin,rBin);\r
253 fHistograms->FillHistogram(nameMCMappingPhiR, daughter0->Vz(), particle->Eta());\r
254 \r
255 TString nameMCMappingPhi="";\r
256 nameMCMappingPhi.Form("MC_EP_Mapping-Phi%02d",phiBin);\r
257 fHistograms->FillHistogram(nameMCMappingPhi, particle->Eta());\r
258 \r
259 TString nameMCMappingR="";\r
260 nameMCMappingR.Form("MC_EP_Mapping-R%02d",rBin);\r
261 fHistograms->FillHistogram(nameMCMappingR, particle->Eta());\r
262 \r
263 TString nameMCMappingPhiInR="";\r
a6f44de5 264 nameMCMappingPhiInR.Form("MC_EP_Mapping_Phi_R-%02d",rBin);\r
265 fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);\r
9fdb2372 266 //end mapping\r
267\r
268 fHistograms->FillHistogram("MC_EP_R",daughter0->R());\r
269 fHistograms->FillHistogram("MC_EP_ZR",daughter0->Vz(),daughter0->R());\r
270 fHistograms->FillHistogram("MC_EP_XY",daughter0->Vx(),daughter0->Vy());\r
271 fHistograms->FillHistogram("MC_EP_OpeningAngle",GetMCOpeningAngle(daughter0, daughter1));\r
a6f44de5 272 }// end if((daughter0->GetPdgCode() == -11 && daughter1->GetPdgCode()) == 11 ||....... approx 20 lines above\r
273 }// end if(particle->GetNDaughters() >= 2){\r
274 } // end if( fStack->Particle(particle->GetMother(0))->GetPdgCode() != 22 )\r
9fdb2372 275 if( fStack->Particle(particle->GetMother(0))->GetPdgCode()==10441 ||//chic0 \r
276 fStack->Particle(particle->GetMother(0))->GetPdgCode()==20443 ||//psi2S\r
277 fStack->Particle(particle->GetMother(0))->GetPdgCode()==445 //chic2\r
a6f44de5 278 ){ \r
9fdb2372 279 fMCGammaChic.push_back(particle);\r
a6f44de5 280 }\r
281 }// end if(particle->GetMother(0)>-1)\r
9fdb2372 282 else{//means we have a primary particle\r
283 fHistograms->FillHistogram("MC_DirectGamma_Energy",particle->Energy());\r
284 fHistograms->FillHistogram("MC_DirectGamma_Pt", particle->Pt());\r
285 fHistograms->FillHistogram("MC_DirectGamma_Eta", particle->Eta());\r
286 fHistograms->FillHistogram("MCDirectGammaPhi", tmpPhi);\r
287\r
288 //adding the conversion points from all gammas with e+e- daughters\r
289 if(particle->GetNDaughters() == 2){\r
290 TParticle* daughter0 = (TParticle*)fStack->Particle(particle->GetFirstDaughter());\r
291 TParticle* daughter1 = (TParticle*)fStack->Particle(particle->GetLastDaughter());\r
292 if((daughter0->GetPdgCode() == -11 && daughter1->GetPdgCode() == 11) ||\r
293 (daughter0->GetPdgCode() == 11 && daughter1->GetPdgCode() == -11)){\r
294 \r
295 fHistograms->FillHistogram("MC_EP_R",daughter0->R());\r
296 fHistograms->FillHistogram("MC_EP_ZR",daughter0->Vz(),daughter0->R());\r
297 fHistograms->FillHistogram("MC_EP_XY",daughter0->Vx(),daughter0->Vy());\r
298 fHistograms->FillHistogram("MC_EP_OpeningAngle",GetMCOpeningAngle(daughter0, daughter1));\r
299\r
300 }\r
301 }\r
a6f44de5 302 }// end else\r
303 }// end if (particle->GetPdgCode()== 22){\r
9fdb2372 304 else if (TMath::Abs(particle->GetPdgCode())== 11){ // Means we have an electron or a positron\r
305 if(particle->GetMother(0)>-1){ // means we have a mother\r
306 if( fStack->Particle(particle->GetMother(0))->GetPdgCode()==22 ){ // Means we have a gamma mother\r
307 if(particle->GetPdgCode() == 11){//electron \r
308 fHistograms->FillHistogram("MC_E_Energy", particle->Energy());\r
309 fHistograms->FillHistogram("MC_E_Pt", particle->Pt());\r
310 fHistograms->FillHistogram("MC_E_Eta", particle->Eta());\r
311 fHistograms->FillHistogram("MC_E_Phi", tmpPhi);\r
312 }\r
313 if(particle->GetPdgCode() == -11){//positron \r
314 fHistograms->FillHistogram("MC_P_Energy", particle->Energy());\r
315 fHistograms->FillHistogram("MC_P_Pt", particle->Pt());\r
316 fHistograms->FillHistogram("MC_P_Eta", particle->Eta());\r
317 fHistograms->FillHistogram("MC_P_Phi", tmpPhi);\r
318 }\r
319 }\r
320 }\r
a6f44de5 321 } // end else if (TMath::Abs(particle->GetPdgCode())== 11)\r
9fdb2372 322 else if(particle->GetNDaughters() == 2){\r
323\r
324 TParticle* daughter0 = (TParticle*)fStack->Particle(particle->GetFirstDaughter());\r
325 TParticle* daughter1 = (TParticle*)fStack->Particle(particle->GetLastDaughter());\r
326 if(daughter0->GetPdgCode() == 22 && daughter1->GetPdgCode() == 22){//check for gamma gamma daughters\r
327 \r
328 if(particle->GetPdgCode()==111){//Pi0\r
329 \r
330 fHistograms->FillHistogram("MC_Pi0_Secondaries_Phi", tmpPhi);\r
331\r
332 if( iTracks >= fStack->GetNprimary()){\r
333 \r
334 fHistograms->FillHistogram("MC_Pi0_Secondaries_Eta", particle->Eta());\r
335\r
336 fHistograms->FillHistogram("MC_Pi0_Secondaries_Phi", tmpPhi);\r
337 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt", particle->Pt());\r
338 fHistograms->FillHistogram("MC_Pi0_Secondaries_Energy", particle->Energy());\r
339 fHistograms->FillHistogram("MC_Pi0_Secondaries_R", particle->R());\r
340 fHistograms->FillHistogram("MC_Pi0_Secondaries_ZR", particle->Vz(),particle->R());\r
a6f44de5 341 fHistograms->FillHistogram("MC_Pi0_Secondaries_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));\r
9fdb2372 342 fHistograms->FillHistogram("MC_Pi0_Secondaries_XY", particle->Vx(),particle->Vy());//only fill from one daughter to avoid multiple filling\r
343 }\r
344 else{\r
345 fHistograms->FillHistogram("MC_Pi0_Eta", particle->Eta());\r
346\r
347 fHistograms->FillHistogram("MC_Pi0_Phi", tmpPhi);\r
348 fHistograms->FillHistogram("MC_Pi0_Pt", particle->Pt());\r
349 fHistograms->FillHistogram("MC_Pi0_Energy", particle->Energy());\r
350 fHistograms->FillHistogram("MC_Pi0_R", particle->R());\r
351 fHistograms->FillHistogram("MC_Pi0_ZR", particle->Vz(),particle->R());\r
a6f44de5 352 fHistograms->FillHistogram("MC_Pi0_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));\r
9fdb2372 353 fHistograms->FillHistogram("MC_Pi0_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling\r
354 }\r
355 }\r
356 else if(particle->GetPdgCode()==221){//Eta\r
357 fHistograms->FillHistogram("MC_Eta_Eta", particle->Eta());\r
358\r
359 fHistograms->FillHistogram("MC_Eta_Phi",tmpPhi);\r
360 fHistograms->FillHistogram("MC_Eta_Pt", particle->Pt());\r
361 fHistograms->FillHistogram("MC_Eta_Energy", particle->Energy());\r
362 fHistograms->FillHistogram("MC_Eta_R", particle->R());\r
363 fHistograms->FillHistogram("MC_Eta_ZR", particle->Vz(),particle->R());\r
a6f44de5 364 fHistograms->FillHistogram("MC_Eta_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));\r
9fdb2372 365 fHistograms->FillHistogram("MC_Eta_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling\r
366 }\r
367 \r
368 //the match data should be filled no matter which mother the gamma-gamma comes from\r
369 fHistograms->FillHistogram("MC_Match_Gamma_R", particle->R());\r
370 fHistograms->FillHistogram("MC_Match_Gamma_ZR", particle->Vz(),particle->R());\r
371 fHistograms->FillHistogram("MC_Match_Gamma_XY", particle->Vx(),particle->Vy());\r
372 fHistograms->FillHistogram("MC_Match_Gamma_Mass", particle->GetCalcMass());\r
373 fHistograms->FillHistogram("MC_Match_Gamma_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));\r
374 fHistograms->FillHistogram("MC_Match_Gamma_Energy", particle->Energy());\r
375 fHistograms->FillHistogram("MC_Match_Gamma_Pt", particle->Pt());\r
376 fHistograms->FillHistogram("MC_Match_Gamma_Eta", particle->Eta());\r
377 fHistograms->FillHistogram("MC_Match_Gamma_Phi",tmpPhi);\r
378 }\r
a6f44de5 379 }// end else if(particle->GetNDaughters() == 2)\r
380 }// end for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++)\r
381} // end ProcessMCData\r
9fdb2372 382\r
383void AliAnalysisTaskGammaConversion::ProcessV0s(){\r
384 // see header file for documentation\r
385\r
386 Int_t nSurvivingV0s=0;\r
387 while(fV0Reader->NextV0()){\r
388 nSurvivingV0s++;\r
389 //-------------------------- filling v0 information -------------------------------------\r
390 fHistograms->FillHistogram("ESD_EP_OpeningAngle", fV0Reader->GetOpeningAngle()); \r
391 fHistograms->FillHistogram("ESD_EP_R", fV0Reader->GetXYRadius());\r
392 fHistograms->FillHistogram("ESD_EP_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());\r
393 fHistograms->FillHistogram("ESD_EP_XY", fV0Reader->GetX(),fV0Reader->GetY());\r
394 \r
395 \r
396 fHistograms->FillHistogram("ESD_E_Energy", fV0Reader->GetNegativeTrackEnergy());\r
397 fHistograms->FillHistogram("ESD_E_Pt", fV0Reader->GetNegativeTrackPt());\r
398 fHistograms->FillHistogram("ESD_E_Eta", fV0Reader->GetNegativeTrackEta());\r
399 fHistograms->FillHistogram("ESD_E_Phi", fV0Reader->GetNegativeTrackPhi());\r
400 \r
401 fHistograms->FillHistogram("ESD_P_Energy", fV0Reader->GetPositiveTrackEnergy());\r
402 fHistograms->FillHistogram("ESD_P_Pt", fV0Reader->GetPositiveTrackPt());\r
403 fHistograms->FillHistogram("ESD_P_Eta", fV0Reader->GetPositiveTrackEta());\r
404 fHistograms->FillHistogram("ESD_P_Phi", fV0Reader->GetPositiveTrackPhi());\r
405 \r
406 fHistograms->FillHistogram("ESD_Gamma_Energy", fV0Reader->GetMotherCandidateEnergy());\r
407 fHistograms->FillHistogram("ESD_Gamma_Pt", fV0Reader->GetMotherCandidatePt());\r
408 fHistograms->FillHistogram("ESD_Gamma_Eta", fV0Reader->GetMotherCandidateEta());\r
409 fHistograms->FillHistogram("ESD_Gamma_Phi", fV0Reader->GetMotherCandidatePhi());\r
410\r
411\r
412 // begin mapping\r
413 Int_t rBin = fHistograms->GetRBin(fV0Reader->GetXYRadius());\r
414 Int_t phiBin = fHistograms->GetPhiBin(fV0Reader->GetNegativeTrackPhi());\r
415 Double_t motherCandidateEta= fV0Reader->GetMotherCandidateEta();\r
416 \r
417 TString nameESDMappingPhiR="";\r
418 nameESDMappingPhiR.Form("ESD_EP_Mapping-Phi%02d-R%02d",phiBin,rBin);\r
419 fHistograms->FillHistogram(nameESDMappingPhiR, fV0Reader->GetZ(), motherCandidateEta);\r
420\r
421 TString nameESDMappingPhi="";\r
422 nameESDMappingPhi.Form("ESD_EP_Mapping-Phi%02d",phiBin);\r
423 fHistograms->FillHistogram(nameESDMappingPhi, fV0Reader->GetZ(), motherCandidateEta);\r
424\r
425 TString nameESDMappingR="";\r
426 nameESDMappingR.Form("ESD_EP_Mapping-R%02d",rBin);\r
427 fHistograms->FillHistogram(nameESDMappingR, fV0Reader->GetZ(), motherCandidateEta); \r
428\r
429 TString nameESDMappingPhiInR="";\r
a6f44de5 430 nameESDMappingPhiInR.Form("ESD_EP_Mapping_Phi_R-%02d",rBin);\r
431 fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());\r
9fdb2372 432 // end mapping\r
433 \r
434 fKFReconstructedGammas.push_back(*fV0Reader->GetMotherCandidateKFCombination());\r
435\r
436 //----------------------------------- checking for "real" conversions (MC match) --------------------------------------\r
437 if(fDoMCTruth){\r
438 if(fV0Reader->HasSameMCMother() == kFALSE){\r
439 continue;\r
440 }\r
a6f44de5 441 \r
9fdb2372 442 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();\r
443 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();\r
a6f44de5 444\r
9fdb2372 445 if(negativeMC->GetPdgCode()!=11 || positiveMC->GetPdgCode()!=-11){\r
446 continue;\r
447 }\r
9fdb2372 448 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){\r
449 fHistograms->FillHistogram("ESD_Match_Gamma_XY", fV0Reader->GetX(),fV0Reader->GetY());\r
450 fHistograms->FillHistogram("ESD_Match_Gamma_OpeningAngle", fV0Reader->GetOpeningAngle());\r
451 fHistograms->FillHistogram("ESD_Match_Gamma_Pt", fV0Reader->GetMotherCandidatePt());\r
452 fHistograms->FillHistogram("ESD_Match_Gamma_Energy", fV0Reader->GetMotherCandidateEnergy());\r
453 fHistograms->FillHistogram("ESD_Match_Gamma_Eta", fV0Reader->GetMotherCandidateEta());\r
454\r
455 fHistograms->FillHistogram("ESD_Match_Gamma_Phi", fV0Reader->GetMotherCandidatePhi());\r
456 fHistograms->FillHistogram("ESD_Match_Gamma_Mass", fV0Reader->GetMotherCandidateMass());\r
457 fHistograms->FillHistogram("ESD_Match_Gamma_Width", fV0Reader->GetMotherCandidateWidth());\r
458 fHistograms->FillHistogram("ESD_Match_Gamma_Chi2", fV0Reader->GetMotherCandidateChi2());\r
459 fHistograms->FillHistogram("ESD_Match_Gamma_NDF", fV0Reader->GetMotherCandidateNDF());\r
460 fHistograms->FillHistogram("ESD_Match_Gamma_R", fV0Reader->GetXYRadius());\r
461 fHistograms->FillHistogram("ESD_Match_Gamma_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());\r
462\r
463 //resolution\r
464 Double_t mcpt = fV0Reader->GetMotherMCParticle()->Pt();\r
465 Double_t esdpt = fV0Reader->GetMotherCandidatePt();\r
466 Double_t resdPt = 0;\r
467 if(mcpt != 0){\r
468 resdPt = ((esdpt - mcpt)/mcpt)*100;\r
469 }\r
470\r
471 fHistograms->FillHistogram("Resolution_dPt", mcpt, resdPt);\r
472 fHistograms->FillHistogram("Resolution_MC_Pt", mcpt);\r
473 fHistograms->FillHistogram("Resolution_ESD_Pt", esdpt);\r
474 \r
475 Double_t resdZ = 0;\r
476 if(fV0Reader->GetNegativeMCParticle()->Vz() != 0){\r
477 resdZ = ((fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz())/fV0Reader->GetNegativeMCParticle()->Vz())*100;\r
478 }\r
479 \r
480 fHistograms->FillHistogram("Resolution_dZ", fV0Reader->GetNegativeMCParticle()->Vz(), resdZ);\r
481 fHistograms->FillHistogram("Resolution_MC_Z", fV0Reader->GetNegativeMCParticle()->Vz());\r
482 fHistograms->FillHistogram("Resolution_ESD_Z", fV0Reader->GetZ());\r
483 \r
484 Double_t resdR = 0;\r
485 if(fV0Reader->GetNegativeMCParticle()->R() != 0){\r
486 resdR = ((fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R())/fV0Reader->GetNegativeMCParticle()->R())*100;\r
487 }\r
a6f44de5 488 fHistograms->FillHistogram("Resolution_dR", fV0Reader->GetNegativeMCParticle()->R(), resdR);\r
9fdb2372 489 fHistograms->FillHistogram("Resolution_MC_R", fV0Reader->GetNegativeMCParticle()->R());\r
490 fHistograms->FillHistogram("Resolution_ESD_R", fV0Reader->GetXYRadius());\r
491 fHistograms->FillHistogram("Resolution_dR_dPt", resdR, resdPt);\r
492 }\r
493 }\r
494 }\r
495 fHistograms->FillHistogram("NumberOfSurvivingV0s", nSurvivingV0s);\r
496 fHistograms->FillHistogram("NumberOfV0s", fV0Reader->GetNumberOfV0s());\r
497}\r
498\r
499void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){\r
500 // see header file for documentation\r
501\r
502 for(UInt_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammas.size();firstGammaIndex++){\r
503 for(UInt_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammas.size();secondGammaIndex++){\r
504 AliKFParticle * twoGammaDecayCandidateDaughter0 = &fKFReconstructedGammas[firstGammaIndex];\r
505 AliKFParticle * twoGammaDecayCandidateDaughter1 = &fKFReconstructedGammas[secondGammaIndex];\r
506\r
507 \r
508 AliKFParticle *twoGammaCandidate = new AliKFParticle(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);\r
509\r
510 Double_t massTwoGammaCandidate =0.;\r
511 Double_t widthTwoGammaCandidate = 0.;\r
512 Double_t chi2TwoGammaCandidate =10000.; \r
513 twoGammaCandidate->GetMass(massTwoGammaCandidate,widthTwoGammaCandidate);\r
514 if(twoGammaCandidate->GetNDF()>0){\r
515 chi2TwoGammaCandidate = twoGammaCandidate->GetChi2()/twoGammaCandidate->GetNDF();\r
a6f44de5 516 if(chi2TwoGammaCandidate>0 && chi2TwoGammaCandidate<fV0Reader->GetChi2CutMeson()){\r
9fdb2372 517\r
518 TVector3 vectorTwoGammaCandidate(twoGammaCandidate->Px(),twoGammaCandidate->Py(),twoGammaCandidate->Pz());\r
519\r
520 Double_t openingAngleTwoGammaCandidate = twoGammaDecayCandidateDaughter0->GetAngle(*twoGammaDecayCandidateDaughter1);\r
521\r
522 //Calculating by hand the radius\r
523 Double_t tmpX= twoGammaCandidate->GetX();\r
524 Double_t tmpY= twoGammaCandidate->GetY();\r
525 \r
526 Double_t radiusTwoGammaCandidate = TMath::Sqrt(tmpX*tmpX + tmpY*tmpY);\r
527\r
a6f44de5 528 fHistograms->FillHistogram("ESD_TwoGammaCombination_GammaDaughter_OpeningAngle", openingAngleTwoGammaCandidate);\r
9fdb2372 529 fHistograms->FillHistogram("ESD_TwoGammaCombination_Energy", twoGammaCandidate->GetE());\r
530 fHistograms->FillHistogram("ESD_TwoGammaCombination_Pt", sqrt(twoGammaCandidate->GetPx()*twoGammaCandidate->GetPx()+twoGammaCandidate->GetPy()*twoGammaCandidate->GetPy()));\r
531 fHistograms->FillHistogram("ESD_TwoGammaCombination_Eta", vectorTwoGammaCandidate.Eta());\r
532 fHistograms->FillHistogram("ESD_TwoGammaCombination_Phi", vectorTwoGammaCandidate.Phi());\r
533 fHistograms->FillHistogram("ESD_TwoGammaCombination_Mass", massTwoGammaCandidate);\r
534 fHistograms->FillHistogram("ESD_TwoGammaCombination_R", radiusTwoGammaCandidate);\r
535 fHistograms->FillHistogram("ESD_TwoGammaCombination_ZR", tmpY, radiusTwoGammaCandidate);\r
a6f44de5 536 fHistograms->FillHistogram("ESD_TwoGammaCombination_XY", tmpX, tmpY);\r
537 fHistograms->FillHistogram("InvMass_vs_Pt_Spectra",massTwoGammaCandidate ,sqrt(twoGammaCandidate->GetPx()*twoGammaCandidate->GetPx()+twoGammaCandidate->GetPy()*twoGammaCandidate->GetPy()));\r
9fdb2372 538 }\r
539 }\r
540 delete twoGammaCandidate;\r
541 }\r
542 }\r
9fdb2372 543}\r
544\r
545void AliAnalysisTaskGammaConversion::CalculateBackground(){\r
546 // see header file for documentation\r
547\r
548 vector<AliKFParticle> vectorCurrentEventGoodV0s = fV0Reader->GetCurrentEventGoodV0s();\r
549 vector<AliKFParticle> vectorPreviousEventGoodV0s = fV0Reader->GetPreviousEventGoodV0s();\r
550 for(UInt_t iCurrent=0;iCurrent<vectorCurrentEventGoodV0s.size();iCurrent++){\r
551 AliKFParticle * currentEventGoodV0 = &vectorCurrentEventGoodV0s.at(iCurrent);\r
552 for(UInt_t iPrevious=0;iPrevious<vectorPreviousEventGoodV0s.size();iPrevious++){\r
553 AliKFParticle * previousGoodV0 = &vectorPreviousEventGoodV0s.at(iPrevious);\r
554\r
555 AliKFParticle *backgroundCandidate = new AliKFParticle(*currentEventGoodV0,*previousGoodV0);\r
556\r
557 Double_t massBG =0.;\r
558 Double_t widthBG = 0.;\r
559 Double_t chi2BG =10000.; \r
560 backgroundCandidate->GetMass(massBG,widthBG);\r
561 if(backgroundCandidate->GetNDF()>0){\r
562 chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();\r
a6f44de5 563 if(chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()){\r
9fdb2372 564\r
565 TVector3 vectorBGCandidate(backgroundCandidate->Px(),backgroundCandidate->Py(),backgroundCandidate->Pz());\r
566\r
567 Double_t openingAngleBG = currentEventGoodV0->GetAngle(*previousGoodV0);\r
568\r
569 //Calculating by hand the radius (find a better way)\r
570 Double_t tmpX= backgroundCandidate->GetX();\r
571 Double_t tmpY= backgroundCandidate->GetY();\r
572 \r
573 Double_t radiusBG = TMath::Sqrt(tmpX*tmpX + tmpY*tmpY);\r
574\r
575 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);\r
576 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());\r
577 fHistograms->FillHistogram("ESD_Background_Pt", sqrt(backgroundCandidate->GetPx()*backgroundCandidate->GetPx()+backgroundCandidate->GetPy()*backgroundCandidate->GetPy()));\r
578 fHistograms->FillHistogram("ESD_Background_Eta", vectorBGCandidate.Eta());\r
579 fHistograms->FillHistogram("ESD_Background_Phi", vectorBGCandidate.Phi());\r
580 fHistograms->FillHistogram("ESD_Background_Mass", massBG);\r
581 fHistograms->FillHistogram("ESD_Background_R", radiusBG);\r
582 fHistograms->FillHistogram("ESD_Background_ZR", tmpY, radiusBG);\r
583 fHistograms->FillHistogram("ESD_Background_XY", tmpX, tmpY);\r
a6f44de5 584 fHistograms->FillHistogram("Background_InvMass_vs_Pt_Spectra",massBG,sqrt(backgroundCandidate->GetPx()*backgroundCandidate->GetPx()+backgroundCandidate->GetPy()*backgroundCandidate->GetPy()));\r
9fdb2372 585 }\r
586 }\r
587 delete backgroundCandidate; \r
588 }\r
589 }\r
590}\r
591\r
592void AliAnalysisTaskGammaConversion::Terminate(Option_t */*option*/)\r
593{\r
594 // Terminate analysis\r
595 //\r
596 AliDebug(1,"Do nothing in Terminate");\r
597}\r
598\r
599void AliAnalysisTaskGammaConversion::UserCreateOutputObjects()\r
600{\r
601 // Create the output container\r
602 if(fOutputContainer != NULL){\r
603 delete fOutputContainer;\r
604 fOutputContainer = NULL;\r
605 }\r
606 if(fOutputContainer == NULL){\r
607 fOutputContainer = new TList();\r
608 }\r
609 fHistograms->GetOutputContainer(fOutputContainer);\r
a6f44de5 610 fOutputContainer->SetName(GetName());\r
9fdb2372 611}\r
612\r
613Double_t AliAnalysisTaskGammaConversion::GetMCOpeningAngle(TParticle* daughter0, TParticle* daughter1) const{\r
614 //helper function\r
615 TVector3 v3D0(daughter0->Px(),daughter0->Py(),daughter0->Pz());\r
616 TVector3 v3D1(daughter1->Px(),daughter1->Py(),daughter1->Pz());\r
617 return v3D0.Angle(v3D1);\r
618}\r