]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/hfe/AliHFEelecbackground.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Added names for the steps
[u/mrichter/AliRoot.git] / PWG3 / hfe / AliHFEelecbackground.cxx
CommitLineData
02524e30 1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15//
16//
17// First implementation of a class
18// to reject tagged electron coming from conversion, pi0 and eta
19// by calculating the e+e- invariant mass
20// of the tagged electron with other tracks
21// after looser cuts for the partner.
22// PostProcess should extract the background yield
23// If running with MC, it can be compared to the expected background yield
24//
25// Authors:
26// Raphaelle Bailhache <rbailhache@ikf.uni-frankfurt.de > <R.Bailhache@gsi.de >
27//
28//
29
30#include <THnSparse.h>
31#include <TParticle.h>
32#include <TFile.h>
33#include <TCanvas.h>
34#include <TAxis.h>
35#include <TH1D.h>
36#include <TString.h>
37#include <TLegend.h>
38
39#include <AliESDEvent.h>
40#include <AliAODEvent.h>
41#include <AliESDtrack.h>
42#include <AliAODTrack.h>
43#include "AliHFEelecbackground.h"
44#include <AliMCEvent.h>
45#include <AliStack.h>
46
47
48ClassImp(AliHFEelecbackground)
49
50 const Double_t AliHFEelecbackground::fgkMe = 0.00051099892;
51
52//___________________________________________________________________________________________
53AliHFEelecbackground::AliHFEelecbackground():
54 fESD1(0x0)
55 ,fAOD1(0x0)
56 ,fMCEvent(0x0)
57 ,fBz(0)
58 ,fkVertex(0x0)
59 ,fUseMCPID(kFALSE)
60 ,fPtESD(0.0)
61 ,fIndexTrack(0)
62 ,fIndexTrackPart(0)
63 ,fPdg(0)
64 ,fLabMother(-1)
65 ,fIsFrom(-1)
66 ,fMotherGamma(-1)
67 ,fMotherPi0(-1)
68 ,fMotherC(-1)
69 ,fMotherB(-1)
70 ,fMotherEta(-1)
71 ,fIsPartner(kFALSE)
72 ,fList(0x0)
73 ,fListPostProcess(0x0)
74{
75 //
76 // Default constructor
77 //
78
79}
80
81//_______________________________________________________________________________________________
82AliHFEelecbackground::AliHFEelecbackground(const AliHFEelecbackground &p):
83 TObject(p)
84 ,fESD1(0x0)
85 ,fAOD1(0x0)
86 ,fMCEvent(0x0)
87 ,fBz(p.fBz)
88 ,fkVertex(p.fkVertex)
89 ,fUseMCPID(p.fUseMCPID)
90 ,fPtESD(p.fPtESD)
91 ,fIndexTrack(0)
92 ,fIndexTrackPart(0)
93 ,fPdg(0)
94 ,fLabMother(-1)
95 ,fIsFrom(-1)
96 ,fMotherGamma(-1)
97 ,fMotherPi0(-1)
98 ,fMotherC(-1)
99 ,fMotherB(-1)
100 ,fMotherEta(-1)
101 ,fIsPartner(kFALSE)
102 ,fList(0x0)
103 ,fListPostProcess(0x0)
104{
105 //
106 // Copy constructor
107 //
108}
109
110//_______________________________________________________________________________________________
111AliHFEelecbackground&
112AliHFEelecbackground::operator=(const AliHFEelecbackground &)
113{
114 //
115 // Assignment operator
116 //
117
118 AliInfo("Not yet implemented.");
119 return *this;
120}
121
122//_______________________________________________________________________________________________
123AliHFEelecbackground::~AliHFEelecbackground()
124{
125 //
126 // Destructor
127 //
128
129 if(fList){
130 fList->Clear();
131 delete fList;
132 }
133
134 if(fListPostProcess){
135 fListPostProcess->Clear();
136 delete fListPostProcess;
137 }
138}
139//___________________________________________________________________________________________
140Bool_t AliHFEelecbackground::Load(const Char_t *filename)
141{
142 //
143 // Generic container loader
144 //
145
146 if(!TFile::Open(filename)){
147 return kFALSE;
148 }
149 TList *o = 0x0;
150 if(!(o = (TList*)gFile->Get("Results"))){
151 return kFALSE;
152 }
153 TList *oe = 0x0;
154 if(!(oe = (TList*)dynamic_cast<TList *>(o->FindObject("HFEelecbackground")))){
155 return kFALSE;
156 }
157 fList = (TList*)oe->Clone("HFEelecbackground");
158 gFile->Close();
159 return kTRUE;
160}
161//_______________________________________________________________________________________________
162void AliHFEelecbackground::CreateHistograms(TList* const qaList)
163{
164 //
165 // create histograms
166 //
167 if(!qaList) return;
168
169 fList = qaList;
170 fList->SetName("HFEelecbackground");
171
172 // bins
173
174 Int_t nBinsPt = 25;
175 Double_t minPt = 0.001;
176 Double_t maxPt = 10.0;
177
178 Int_t nBinsInv = 50;
179 Double_t minInv = 0.0;
180 Double_t maxInv = 0.2;
181
182 Int_t nBinsOp = 50;
183 Double_t minOp = 0.0;
184 Double_t maxOp = 2;
185
186 Double_t *binLimLogPt = new Double_t[nBinsPt+1];
187 Double_t *binLimPt = new Double_t[nBinsPt+1];
188 for(Int_t i=0; i<=nBinsPt; i++) binLimLogPt[i]=(Double_t)TMath::Log10(minPt) + (TMath::Log10(maxPt)-TMath::Log10(minPt))/nBinsPt*(Double_t)i ;
189 for(Int_t i=0; i<=nBinsPt; i++) binLimPt[i]=(Double_t)TMath::Power(10,binLimLogPt[i]);
190
191 Double_t *binLimInv = new Double_t[nBinsInv+1];
192 for(Int_t i=0; i<=nBinsInv; i++) binLimInv[i]=(Double_t)minInv + (maxInv-minInv) /nBinsInv*(Double_t)i ;
193
194 Double_t *binLimOp = new Double_t[nBinsOp+1];
195 for(Int_t i=0; i<=nBinsOp; i++) binLimOp[i]=(Double_t)minOp + (maxOp-minOp) /nBinsOp*(Double_t)i ;
196
197
198 const Int_t nvarO = 3; // ptrectaggede, ptrecmother, openingangle or invmass
199
200 Int_t iBinOInv[nvarO];
201 iBinOInv[0]=nBinsPt;
202 iBinOInv[1]=nBinsPt;
203 iBinOInv[2]=nBinsInv;
204
205 Int_t iBinOOp[nvarO];
206 iBinOOp[0]=nBinsPt;
207 iBinOOp[1]=nBinsPt;
208 iBinOOp[2]=nBinsOp;
209
210 //
211 //
212 //
213
214 THnSparseF *openinganglepp = new THnSparseF("openinganglepp","",nvarO,iBinOOp);
215 openinganglepp->SetBinEdges(0,binLimPt);
216 openinganglepp->SetBinEdges(1,binLimPt);
217 openinganglepp->SetBinEdges(2,binLimOp);
218 openinganglepp->Sumw2();
219
220 THnSparseF *openinganglenn = new THnSparseF("openinganglenn","",nvarO,iBinOOp);
221 openinganglenn->SetBinEdges(0,binLimPt);
222 openinganglenn->SetBinEdges(1,binLimPt);
223 openinganglenn->SetBinEdges(2,binLimOp);
224 openinganglenn->Sumw2();
225
226 THnSparseF *openingangless = new THnSparseF("openingangless","",nvarO,iBinOOp);
227 openingangless->SetBinEdges(0,binLimPt);
228 openingangless->SetBinEdges(1,binLimPt);
229 openingangless->SetBinEdges(2,binLimOp);
230 openingangless->Sumw2();
231
232 THnSparseF *openingangler = new THnSparseF("openingangler","",nvarO,iBinOOp);
233 openingangler->SetBinEdges(0,binLimPt);
234 openingangler->SetBinEdges(1,binLimPt);
235 openingangler->SetBinEdges(2,binLimOp);
236 openingangler->Sumw2();
237
238 THnSparseF *openingangleos = new THnSparseF("openingangleos","",nvarO,iBinOOp);
239 openingangleos->SetBinEdges(0,binLimPt);
240 openingangleos->SetBinEdges(1,binLimPt);
241 openingangleos->SetBinEdges(2,binLimOp);
242 openingangleos->Sumw2();
243
244 THnSparseF *openinganglepi0=0x0;
245 THnSparseF *openingangleeta=0x0;
246 THnSparseF *openinganglegamma=0x0;
247 THnSparseF *openingangleC=0x0;
248 THnSparseF *openingangleB=0x0;
249
250 if(HasMCData()) {
251
252 openinganglepi0 = new THnSparseF("openinganglepi0","",nvarO,iBinOOp);
253 openinganglepi0->SetBinEdges(0,binLimPt);
254 openinganglepi0->SetBinEdges(1,binLimPt);
255 openinganglepi0->SetBinEdges(2,binLimOp);
256 openinganglepi0->Sumw2();
257
258 openingangleeta = new THnSparseF("openingangleeta","",nvarO,iBinOOp);
259 openingangleeta->SetBinEdges(0,binLimPt);
260 openingangleeta->SetBinEdges(1,binLimPt);
261 openingangleeta->SetBinEdges(2,binLimOp);
262 openingangleeta->Sumw2();
263
264 openinganglegamma = new THnSparseF("openinganglegamma","",nvarO,iBinOOp);
265 openinganglegamma->SetBinEdges(0,binLimPt);
266 openinganglegamma->SetBinEdges(1,binLimPt);
267 openinganglegamma->SetBinEdges(2,binLimOp);
268 openinganglegamma->Sumw2();
269
270 openingangleC = new THnSparseF("openingangleC","",nvarO,iBinOOp);
271 openingangleC->SetBinEdges(0,binLimPt);
272 openingangleC->SetBinEdges(1,binLimPt);
273 openingangleC->SetBinEdges(2,binLimOp);
274 openingangleC->Sumw2();
275
276 openingangleB = new THnSparseF("openingangleB","",nvarO,iBinOOp);
277 openingangleB->SetBinEdges(0,binLimPt);
278 openingangleB->SetBinEdges(1,binLimPt);
279 openingangleB->SetBinEdges(2,binLimOp);
280 openingangleB->Sumw2();
281
282 }
283
284 //
285
286 THnSparseF *invmasspp = new THnSparseF("invmasspp","",nvarO,iBinOInv);
287 invmasspp->SetBinEdges(0,binLimPt);
288 invmasspp->SetBinEdges(1,binLimPt);
289 invmasspp->SetBinEdges(2,binLimInv);
290 invmasspp->Sumw2();
291
292 THnSparseF *invmassnn = new THnSparseF("invmassnn","",nvarO,iBinOInv);
293 invmassnn->SetBinEdges(0,binLimPt);
294 invmassnn->SetBinEdges(1,binLimPt);
295 invmassnn->SetBinEdges(2,binLimInv);
296 invmassnn->Sumw2();
297
298 THnSparseF *invmassss = new THnSparseF("invmassss","",nvarO,iBinOInv);
299 invmassss->SetBinEdges(0,binLimPt);
300 invmassss->SetBinEdges(1,binLimPt);
301 invmassss->SetBinEdges(2,binLimInv);
302 invmassss->Sumw2();
303
304 THnSparseF *invmassr = new THnSparseF("invmassr","",nvarO,iBinOInv);
305 invmassr->SetBinEdges(0,binLimPt);
306 invmassr->SetBinEdges(1,binLimPt);
307 invmassr->SetBinEdges(2,binLimInv);
308 invmassr->Sumw2();
309
310 THnSparseF *invmassos = new THnSparseF("invmassos","",nvarO,iBinOInv);
311 invmassos->SetBinEdges(0,binLimPt);
312 invmassos->SetBinEdges(1,binLimPt);
313 invmassos->SetBinEdges(2,binLimInv);
314 invmassos->Sumw2();
315
316 THnSparseF *invmasspi0=0x0;
317 THnSparseF *invmasseta=0x0;
318 THnSparseF *invmassgamma=0x0;
319 THnSparseF *invmassC=0x0;
320 THnSparseF *invmassB=0x0;
321
322 if(HasMCData()) {
323
324 invmasspi0 = new THnSparseF("invmasspi0","",nvarO,iBinOInv);
325 invmasspi0->SetBinEdges(0,binLimPt);
326 invmasspi0->SetBinEdges(1,binLimPt);
327 invmasspi0->SetBinEdges(2,binLimInv);
328 invmasspi0->Sumw2();
329
330 invmasseta = new THnSparseF("invmasseta","",nvarO,iBinOInv);
331 invmasseta->SetBinEdges(0,binLimPt);
332 invmasseta->SetBinEdges(1,binLimPt);
333 invmasseta->SetBinEdges(2,binLimInv);
334 invmasseta->Sumw2();
335
336 invmassgamma = new THnSparseF("invmassgamma","",nvarO,iBinOInv);
337 invmassgamma->SetBinEdges(0,binLimPt);
338 invmassgamma->SetBinEdges(1,binLimPt);
339 invmassgamma->SetBinEdges(2,binLimInv);
340 invmassgamma->Sumw2();
341
342 invmassC = new THnSparseF("invmassC","",nvarO,iBinOInv);
343 invmassC->SetBinEdges(0,binLimPt);
344 invmassC->SetBinEdges(1,binLimPt);
345 invmassC->SetBinEdges(2,binLimInv);
346 invmassC->Sumw2();
347
348 invmassB = new THnSparseF("invmassB","",nvarO,iBinOInv);
349 invmassB->SetBinEdges(0,binLimPt);
350 invmassB->SetBinEdges(1,binLimPt);
351 invmassB->SetBinEdges(2,binLimInv);
352 invmassB->Sumw2();
353
354 }
355
356 //
357 //
358 //
359
360 fList->AddAt(openinganglepp,kPp);
361 fList->AddAt(openinganglenn,kNn);
362 fList->AddAt(openingangless,kSs);
363 fList->AddAt(openingangler,kR);
364 fList->AddAt(openingangleos,kOs);
365 if(HasMCData()) {
366 fList->AddAt(openinganglepi0,2*kNSignComb+kElectronFromPi0);
367 fList->AddAt(openingangleeta,2*kNSignComb+kElectronFromEta);
368 fList->AddAt(openinganglegamma,2*kNSignComb+kElectronFromGamma);
369 fList->AddAt(openingangleC,2*kNSignComb+kElectronFromC);
370 fList->AddAt(openingangleB,2*kNSignComb+kElectronFromB);
371 }
372
373 fList->AddAt(invmasspp,kNSignComb+kPp);
374 fList->AddAt(invmassnn,kNSignComb+kNn);
375 fList->AddAt(invmassss,kNSignComb+kSs);
376 fList->AddAt(invmassr,kNSignComb+kR);
377 fList->AddAt(invmassos,kNSignComb+kOs);
378 if(HasMCData()) {
379 fList->AddAt(invmasspi0,2*kNSignComb+kNMCInfo+kElectronFromPi0);
380 fList->AddAt(invmasseta,2*kNSignComb+kNMCInfo+kElectronFromEta);
381 fList->AddAt(invmassgamma,2*kNSignComb+kNMCInfo+kElectronFromGamma);
382 fList->AddAt(invmassC,2*kNSignComb+kNMCInfo+kElectronFromC);
383 fList->AddAt(invmassB,2*kNSignComb+kNMCInfo+kElectronFromB);
384 }
385
386}
387//_______________________________________________________________________________________________
388void AliHFEelecbackground::PairAnalysis(AliESDtrack* const track, AliESDtrack* const trackPart)
389{
390 //
391 // calculate (tagged e-partner) dca, opening angle, invariant mass
392 //
393
394 ////////////////////////
395 // Partner track cut
396 ////////////////////////
397 if(!SingleTrackCut(trackPart)) return;
398
399 ////////////////////////
400 // Take label
401 ////////////////////////
402 Int_t indexTrack = TMath::Abs(track->GetLabel());
403 Int_t indexTrackPart = TMath::Abs(trackPart->GetLabel());
404
405 /////////////////////////
406 // If MC data
407 ////////////////////////
408
409 if(HasMCData()) {
410
411 // Take info track if not already done
412 if(indexTrack!= fIndexTrack) {
413 fIsFrom = -1;
414
415 fPdg = GetPdg(indexTrack);
416 fLabMother = GetLabMother(indexTrack);
417
418 fMotherGamma = IsMotherGamma(indexTrack);
419 if((fMotherGamma != -1) && ((TMath::Abs(fPdg)) == 11)) fIsFrom = kElectronFromGamma;
420 fMotherPi0 = IsMotherPi0(indexTrack);
421 if((fMotherPi0 != -1) && ((TMath::Abs(fPdg)) == 11)) fIsFrom = kElectronFromPi0;
422 fMotherC = IsMotherC(indexTrack);
423 if((fMotherC != -1) && ((TMath::Abs(fPdg)) == 11)) fIsFrom = kElectronFromC;
424 fMotherB = IsMotherB(indexTrack);
425 if((fMotherB != -1) && ((TMath::Abs(fPdg)) == 11)) fIsFrom = kElectronFromB;
426 fMotherEta = IsMotherEta(indexTrack);
427 if((fMotherEta != -1) && ((TMath::Abs(fPdg)) == 11)) fIsFrom = kElectronFromEta;
428
429 }
430
431 // Look at trackPart
432 Int_t pdgPart = GetPdg(indexTrackPart);
433 if(TMath::Abs(pdgPart) == 11) {
434 Int_t labMotherPart = GetLabMother(indexTrackPart);
435 if((labMotherPart == fLabMother) && (indexTrack != indexTrackPart) && (TMath::Abs(pdgPart)==11) && (fPdg*pdgPart < 0) && (fLabMother >=0) && (fLabMother < (((AliStack *)fMCEvent->Stack())->GetNtrack()))) fIsPartner = kTRUE;
436 }
437
438 }
439
440 //////////////////////
441 // Sign
442 /////////////////////
443 Int_t sign = -1;
444 if((track->Charge() > 0.0) && (trackPart->Charge() > 0.0)) sign = kPp;
445 if((track->Charge() < 0.0) && (trackPart->Charge() < 0.0)) sign = kNn;
446 if(((track->Charge() > 0.0) && (trackPart->Charge() < 0.0)) || ((track->Charge() < 0.0) && (trackPart->Charge() > 0.0))) sign = kOs;
447
448 //////////////////////
449 // DCA
450 /////////////////////
451
452 Double_t xthis,xp;
453 Double_t dca = track->GetDCA(trackPart,fBz,xthis,xp);
454 if(TMath::Abs(dca) > 3.0) return;
455
456 /////////////////////////////
457 // Propagate
458 ////////////////////////////
459
460 Double_t norradius = TMath::Sqrt(fkVertex->GetX()*fkVertex->GetX()+fkVertex->GetY()*fkVertex->GetY());
461
462 AliESDtrack *trackCopy = new AliESDtrack(*track);
463 AliESDtrack *trackPartCopy = new AliESDtrack(*trackPart);
464 Bool_t propagateok = kTRUE;
465 if((!(trackPartCopy->PropagateTo(norradius,fBz))) || (!(trackCopy->PropagateTo(norradius,fBz)))) propagateok = kFALSE;
466 if(!propagateok) {
467 if(trackCopy) delete trackCopy;
468 if(trackPartCopy) delete trackPartCopy;
469 return;
470 }
471
472 ///////////////////////////////////
473 // Calcul mother variables
474 ///////////////////////////////////
475 Double_t results[5];
476 Double_t resultsr[5];
477
478 CalculateMotherVariable(trackCopy,trackPartCopy,&results[0]);
479 CalculateMotherVariableR(trackCopy,trackPartCopy,&resultsr[0]);
480
481 /////////////////////////////////////
482 // Fill
483 /////////////////////////////////////
484
485 FillOutput(results, resultsr, sign);
486
487 if(trackCopy) delete trackCopy;
488 if(trackPartCopy) delete trackPartCopy;
489
490}
491//_____________________________________________________________________________________
492void AliHFEelecbackground::CalculateMotherVariable(AliESDtrack* const track, AliESDtrack* const trackpart, Double_t *results)
493{
494 //
495 // variables mother and take the pt of the track
496 //
497 // results contain: ptmother, invmass, etamother, phimother, openingangle
498 //
499
500 TVector3 v3Dtagged;
501 TVector3 v3Dpart;
502
503 Double_t *pxyz = new Double_t[3];
504 track->PxPyPz(pxyz);
505 v3Dtagged.SetXYZ(pxyz[0],pxyz[1],pxyz[2]);
506 fPtESD = TMath::Sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]);
507
508 Double_t *pxyzpart = new Double_t[3];
509 trackpart->PxPyPz(pxyzpart);
510 v3Dpart.SetXYZ(pxyzpart[0],pxyzpart[1],pxyzpart[2]);
511
512
513
514 TVector3 motherrec = v3Dtagged + v3Dpart;
515
516 Double_t etaESDmother = motherrec.Eta();
517 Double_t ptESDmother = motherrec.Pt();
518 Double_t phiESDmother = motherrec.Phi();
519 if(phiESDmother > TMath::Pi()) phiESDmother = phiESDmother - (2*TMath::Pi());
520
521
522 // openinganglepropagated
523 Double_t openingangle = v3Dtagged.Angle(v3Dpart);
524
525 // invmass
526 Double_t pESD = TMath::Sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]+pxyz[2]*pxyz[2]);
527 Double_t pESDpart = TMath::Sqrt(pxyzpart[0]*pxyzpart[0]+pxyzpart[1]*pxyzpart[1]+pxyzpart[2]*pxyzpart[2]);
528
529 // e propagate
530 Double_t eESD = TMath::Sqrt(pESD*pESD+fgkMe*fgkMe);
531 Double_t eESDpart = TMath::Sqrt(pESDpart*pESDpart+fgkMe*fgkMe);
532
533 Double_t invmass = TMath::Sqrt((eESD+eESDpart)*(eESD+eESDpart)-(motherrec.Px()*motherrec.Px()+motherrec.Py()*motherrec.Py()+motherrec.Pz()*motherrec.Pz()));
534
535 if(!results) {
536 results = new Double_t[5];
537 }
538 results[0] = ptESDmother;
539 results[1] = etaESDmother;
540 results[2] = phiESDmother;
541 results[3] = invmass;
542 results[4] = openingangle;
543
544}
545//_____________________________________________________________________________________
546void AliHFEelecbackground::CalculateMotherVariableR(AliESDtrack* const track, AliESDtrack* const trackpart, Double_t *results)
547{
548 //
549 // variables mother
550 //
551 // results contain: ptmother, invmass, etamother, phimother, openingangle
552 //
553
554 TVector3 v3Dtagged;
555 TVector3 v3Dpart;
556
557 Double_t *pxyz = new Double_t[3];
558 track->PxPyPz(pxyz);
559 v3Dtagged.SetXYZ(pxyz[0],pxyz[1],pxyz[2]);
560 Double_t *pxyzpart = new Double_t[3];
561 trackpart->PxPyPz(pxyzpart);
562 v3Dpart.SetXYZ(pxyzpart[0],pxyzpart[1],pxyzpart[2]);
563
564 // rotate the partner
565 v3Dpart.RotateZ(TMath::Pi());
566 v3Dpart.GetXYZ(pxyzpart);
567
568
569 TVector3 motherrec = v3Dtagged + v3Dpart;
570
571 Double_t etaESDmother = motherrec.Eta();
572 Double_t ptESDmother = motherrec.Pt();
573 Double_t phiESDmother = motherrec.Phi();
574 if(phiESDmother > TMath::Pi()) phiESDmother = phiESDmother - (2*TMath::Pi());
575
576
577 // openinganglepropagated
578 Double_t openingangle = v3Dtagged.Angle(v3Dpart);
579 //printf("Openingangle %f\n",openingangle);
580
581 // invmass
582 Double_t pESD = TMath::Sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]+pxyz[2]*pxyz[2]);
583 Double_t pESDpart = TMath::Sqrt(pxyzpart[0]*pxyzpart[0]+pxyzpart[1]*pxyzpart[1]+pxyzpart[2]*pxyzpart[2]);
584 // e propagate
585 Double_t eESD = TMath::Sqrt(pESD*pESD+fgkMe*fgkMe);
586 Double_t eESDpart = TMath::Sqrt(pESDpart*pESDpart+fgkMe*fgkMe);
587
588 Double_t invmass = TMath::Sqrt((eESD+eESDpart)*(eESD+eESDpart)-(motherrec.Px()*motherrec.Px()+motherrec.Py()*motherrec.Py()+motherrec.Pz()*motherrec.Pz()));
589
590 if(!results) {
591 results = new Double_t[5];
592 }
593 results[0] = ptESDmother;
594 results[1] = etaESDmother;
595 results[2] = phiESDmother;
596 results[3] = invmass;
597 results[4] = openingangle;
598
599
600}
601//_________________________________________________________________________________
602void AliHFEelecbackground::FillOutput(Double_t *results, Double_t *resultsr, Int_t sign)
603{
604 //
605 // Fill the invariant mass and opening angle distributions
606 //
607
608 Double_t co[3];
609 co[0] = fPtESD;
610
611 switch(sign){
612
613 case kPp:
614 co[1] = results[0];
615 co[2] = TMath::Abs(results[4]);
616 (dynamic_cast<THnSparseF *>(fList->At(kPp)))->Fill(co);
617 (dynamic_cast<THnSparseF *>(fList->At(kSs)))->Fill(co);
618
619
620 co[2] = results[3];
621 if(TMath::Abs(results[4]) < 0.8){
622 (dynamic_cast<THnSparseF *>(fList->At(kPp+kNSignComb)))->Fill(co);
623 (dynamic_cast<THnSparseF *>(fList->At(kSs+kNSignComb)))->Fill(co);
624 }
625
626 break;
627
628 case kNn:
629 co[1] = results[0];
630 co[2] = TMath::Abs(results[4]);
631 (dynamic_cast<THnSparseF *>(fList->At(kNn)))->Fill(co);
632 (dynamic_cast<THnSparseF *>(fList->At(kSs)))->Fill(co);
633
634 co[2] = results[3];
635 if(TMath::Abs(results[4]) < 0.8){
636 (dynamic_cast<THnSparseF *>(fList->At(kNn+kNSignComb)))->Fill(co);
637 (dynamic_cast<THnSparseF *>(fList->At(kSs+kNSignComb)))->Fill(co);
638 }
639 break;
640
641 case kOs:
642 co[1] = results[0];
643 co[2] = TMath::Abs(results[4]);
644 (dynamic_cast<THnSparseF *>(fList->At(kOs)))->Fill(co);
645 if(HasMCData()) {
646 if((fIsFrom == kElectronFromPi0) && (fIsPartner)) (dynamic_cast<THnSparseF *>(fList->At(2*kNSignComb+kElectronFromPi0)))->Fill(co);
647 if((fIsFrom == kElectronFromEta) && (fIsPartner)) (dynamic_cast<THnSparseF *>(fList->At(2*kNSignComb+kElectronFromEta)))->Fill(co);
648 if((fIsFrom == kElectronFromGamma) && (fIsPartner)) (dynamic_cast<THnSparseF *>(fList->At(2*kNSignComb+kElectronFromGamma)))->Fill(co);
649 if((fIsFrom == kElectronFromC) && (fIsPartner)) (dynamic_cast<THnSparseF *>(fList->At(2*kNSignComb+kElectronFromC)))->Fill(co);
650 if((fIsFrom == kElectronFromB) && (fIsPartner)) (dynamic_cast<THnSparseF *>(fList->At(2*kNSignComb+kElectronFromB)))->Fill(co);
651 }
652
653 co[2] = results[3];
654 if(TMath::Abs(results[4]) < 0.8){
655 (dynamic_cast<THnSparseF *>(fList->At(kOs+kNSignComb)))->Fill(co);
656 if(HasMCData()) {
657 if((fIsFrom == kElectronFromPi0) && (fIsPartner)) (dynamic_cast<THnSparseF *>(fList->At(2*kNSignComb+kElectronFromPi0+kNMCInfo)))->Fill(co);
658 if((fIsFrom == kElectronFromEta) && (fIsPartner)) (dynamic_cast<THnSparseF *>(fList->At(2*kNSignComb+kElectronFromEta+kNMCInfo)))->Fill(co);
659 if((fIsFrom == kElectronFromGamma) && (fIsPartner)) (dynamic_cast<THnSparseF *>(fList->At(2*kNSignComb+kElectronFromGamma+kNMCInfo)))->Fill(co);
660 if((fIsFrom == kElectronFromC) && (fIsPartner)) (dynamic_cast<THnSparseF *>(fList->At(2*kNSignComb+kElectronFromC+kNMCInfo)))->Fill(co);
661 if((fIsFrom == kElectronFromB) && (fIsPartner)) (dynamic_cast<THnSparseF *>(fList->At(2*kNSignComb+kElectronFromB+kNMCInfo)))->Fill(co);
662 }
663 }
664
665 // rotated
666 co[1] = resultsr[0];
667 co[2] = TMath::Abs(resultsr[4]);
668
669 (dynamic_cast<THnSparseF *>(fList->At(kR)))->Fill(co);
670
671 co[2] = resultsr[3];
672 if(TMath::Abs(resultsr[4]) < 0.8){
673 (dynamic_cast<THnSparseF *>(fList->At(kR+kNSignComb)))->Fill(co);
674 }
675 break;
676
677 default:
678
679 break;
680
681 }
682}
683//_______________________________________________________________________________________________
684Bool_t AliHFEelecbackground::SingleTrackCut(AliESDtrack* const track) const
685{
686 //
687 // Return minimum quality for the partner
688 //
689
690 UInt_t status = track->GetStatus();
691
692 if(((status & AliESDtrack::kTPCin)==0) && (status & AliESDtrack::kITSin)) {
693
694 Int_t nbcl = track->GetITSclusters(0);
695 if(nbcl > 1) return kTRUE;
696 else return kFALSE;
697
698 }
699
700 if(status & AliESDtrack::kTPCin) {
701
702 if(status & AliESDtrack::kTPCrefit) return kTRUE;
703 else return kFALSE;
704
705 }
706
707 return kFALSE;
708
709}
710//______________________________________________________________________________________________
711void AliHFEelecbackground::SetEvent(AliESDEvent* const ESD)
712{
713 //
714 // Set the AliESD Event, the magnetic field and the primary vertex
715 //
716
717 fESD1=ESD;
718 fBz=fESD1->GetMagneticField();
719 fkVertex = fESD1->GetPrimaryVertex();
720
721}
722//____________________________________________________________________________________________________________
723Int_t AliHFEelecbackground::IsMotherGamma(Int_t tr) {
724
725 //
726 // Return the lab of gamma mother or -1 if not gamma
727 //
728
729 AliStack* stack = fMCEvent->Stack();
730 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
731
732 // Take mother
733 TParticle * particle = stack->Particle(tr);
734 if(!particle) return -1;
735 Int_t imother = particle->GetFirstMother();
736 if((imother < 0) || (imother >= stack->GetNtrack())) return -1;
737 TParticle * mother = stack->Particle(imother);
738 if(!mother) return -1;
739
740 // Check gamma
741 Int_t pdg = mother->GetPdgCode();
742 if(TMath::Abs(pdg) == 22) return imother;
743 if(TMath::Abs(pdg) == 11) {
744 return IsMotherGamma(imother);
745 }
746 return -1;
747
748}
749//
750//____________________________________________________________________________________________________________
751Int_t AliHFEelecbackground::IsMotherPi0(Int_t tr) {
752
753 //
754 // Return the lab of pi0 mother or -1 if not pi0
755 //
756
757 AliStack* stack = fMCEvent->Stack();
758 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
759
760 // Take mother
761 TParticle * particle = stack->Particle(tr);
762 if(!particle) return -1;
763 Int_t imother = particle->GetFirstMother();
764 if((imother < 0) || (imother >= stack->GetNtrack())) return -1;
765 TParticle * mother = stack->Particle(imother);
766 if(!mother) return -1;
767
768 // Check gamma
769 Int_t pdg = mother->GetPdgCode();
770 if(TMath::Abs(pdg) == 111) return imother;
771 if(TMath::Abs(pdg) == 11) {
772 return IsMotherPi0(imother);
773 }
774 return -1;
775
776}
777//____________________________________________________________________________________________________________
778Int_t AliHFEelecbackground::IsMotherEta(Int_t tr) {
779
780 //
781 // Return the lab of pi0 mother or -1 if not pi0
782 //
783
784 AliStack* stack = fMCEvent->Stack();
785 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
786
787 // Take mother
788 TParticle * particle = stack->Particle(tr);
789 if(!particle) return -1;
790 Int_t imother = particle->GetFirstMother();
791 if((imother < 0) || (imother >= stack->GetNtrack())) return -1;
792 TParticle * mother = stack->Particle(imother);
793 if(!mother) return -1;
794
795 // Check gamma
796 Int_t pdg = mother->GetPdgCode();
797 if(TMath::Abs(pdg) == 221) return imother;
798 if(TMath::Abs(pdg) == 11) {
799 return IsMotherEta(imother);
800 }
801 return -1;
802
803}
804//____________________________________________________________________________________________________________
805Int_t AliHFEelecbackground::IsMotherC(Int_t tr) {
806
807 //
808 // Return the lab of signal mother or -1 if not signal
809 //
810
811 AliStack* stack = fMCEvent->Stack();
812 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
813
814 // Take mother
815 TParticle * particle = stack->Particle(tr);
816 if(!particle) return -1;
817 Int_t imother = particle->GetFirstMother();
818 if((imother < 0) || (imother >= stack->GetNtrack())) return -1;
819 TParticle * mother = stack->Particle(imother);
820 if(!mother) return -1;
821
822 // Check gamma
823 Int_t pdg = mother->GetPdgCode();
824 if((TMath::Abs(pdg)==411) || (TMath::Abs(pdg)==421) || (TMath::Abs(pdg)==431) || (TMath::Abs(pdg)==4122) || (TMath::Abs(pdg)==4132) || (TMath::Abs(pdg)==4232) || (TMath::Abs(pdg)==43320)) return imother;
825 if(TMath::Abs(pdg) == 11) {
826 return IsMotherC(imother);
827 }
828 return -1;
829
830}
831//____________________________________________________________________________________________________________
832Int_t AliHFEelecbackground::IsMotherB(Int_t tr) {
833
834 //
835 // Return the lab of signal mother or -1 if not signal
836 //
837
838 AliStack* stack = fMCEvent->Stack();
839 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
840
841 // Take mother
842 TParticle * particle = stack->Particle(tr);
843 if(!particle) return -1;
844 Int_t imother = particle->GetFirstMother();
845 if((imother < 0) || (imother >= stack->GetNtrack())) return -1;
846 TParticle * mother = stack->Particle(imother);
847 if(!mother) return -1;
848
849 // Check gamma
850 Int_t pdg = mother->GetPdgCode();
851 if((TMath::Abs(pdg)==511) || (TMath::Abs(pdg)==521) || (TMath::Abs(pdg)==531) || (TMath::Abs(pdg)==5122) || (TMath::Abs(pdg)==5132) || (TMath::Abs(pdg)==5232) || (TMath::Abs(pdg)==53320)) return imother;
852 if(TMath::Abs(pdg) == 11) {
853 return IsMotherB(imother);
854 }
855 return -1;
856
857}
858//____________________________________________________________________________________________________________
859Int_t AliHFEelecbackground::GetPdg(Int_t tr) {
860
861 //
862 // Simply pdg code
863 //
864
865 AliStack* stack = fMCEvent->Stack();
866 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
867
868 // MC Information
869 TParticle * particle = stack->Particle(tr);
870 if(!particle) return -1;
871 Int_t pdg = particle->GetPdgCode();
872
873 return pdg;
874
875}
876//____________________________________________________________________________________________________________
877Int_t AliHFEelecbackground::GetLabMother(Int_t tr) {
878
879 //
880 // Simply lab mother
881 //
882
883 AliStack* stack = fMCEvent->Stack();
884 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
885
886 // MC Information
887 TParticle * particle = stack->Particle(tr);
888 if(!particle) return -1;
889 Int_t imother = particle->GetFirstMother();
890
891 return imother;
892
893}
894//_______________________________________________________________________________________________
895void AliHFEelecbackground::PostProcess()
896{
897 //
898 // Post process the histos and extract the background pt spectra
899 //
900
901 if(!fList) return;
902
903 // invariant mass input spectra
904 THnSparseF *invmassss = dynamic_cast<THnSparseF *>(fList->FindObject("invmassss"));
905 THnSparseF *invmassr = dynamic_cast<THnSparseF *>(fList->FindObject("invmassr"));
906 THnSparseF *invmassos = dynamic_cast<THnSparseF *>(fList->FindObject("invmassos"));
907 THnSparseF *invmassgamma = dynamic_cast<THnSparseF *>(fList->FindObject("invmassgamma"));
908 THnSparseF *invmasspi0 = dynamic_cast<THnSparseF *>(fList->FindObject("invmasspi0"));
909 THnSparseF *invmasseta = dynamic_cast<THnSparseF *>(fList->FindObject("invmasseta"));
910 THnSparseF *invmassC = dynamic_cast<THnSparseF *>(fList->FindObject("invmassC"));
911 THnSparseF *invmassB = dynamic_cast<THnSparseF *>(fList->FindObject("invmassB"));
912
913 TAxis *ptaxisinvmass = invmassss->GetAxis(0);
914 Int_t nbinsptinvmass = ptaxisinvmass->GetNbins();
915
916 // outputs
917 TH1D **invmassosptproj = new TH1D*[nbinsptinvmass];
918 TH1D **invmassssptproj = new TH1D*[nbinsptinvmass];
919 TH1D **invmassrptproj = new TH1D*[nbinsptinvmass];
920 TH1D **invmassdiffptproj = new TH1D*[nbinsptinvmass];
921 TH1D **invmassgammaptproj = new TH1D*[nbinsptinvmass];
922 TH1D **invmasspi0ptproj = new TH1D*[nbinsptinvmass];
923 TH1D **invmassetaptproj = new TH1D*[nbinsptinvmass];
924 TH1D **invmassCptproj = new TH1D*[nbinsptinvmass];
925 TH1D **invmassBptproj = new TH1D*[nbinsptinvmass];
926
927 TH1D *yieldPtFound = (TH1D *) invmassss->Projection(0);
928 yieldPtFound->SetName("Found yield");
929 yieldPtFound->Reset();
930
931 TH1D *yieldPtSourcesMC = 0x0;
932 if(invmasspi0 && invmasseta && invmassgamma) {
933 yieldPtSourcesMC = (TH1D *) invmassss->Projection(0);
934 yieldPtSourcesMC->SetName("Found yield");
935 yieldPtSourcesMC->Reset();
936 }
937
938 TH1D *yieldPtSignalCutMC = 0x0;
939 if(invmassC && invmassB) {
940 yieldPtSignalCutMC = (TH1D *) invmassss->Projection(0);
941 yieldPtSignalCutMC->SetName("Found yield");
942 yieldPtSignalCutMC->Reset();
943 }
944
945 // canvas
946 Int_t nbrow = (Int_t) (nbinsptinvmass/5);
947 TString namecanvas("Invmassnamecanvas");
948 TCanvas * canvas =new TCanvas(namecanvas,namecanvas,800,800);
949 canvas->Divide(5,nbrow+1);
950
951
952 // loop on pt bins
953 for(Int_t k=1; k <= nbinsptinvmass; k++){
954
955 Double_t lowedge = ptaxisinvmass->GetBinLowEdge(k);
956 Double_t upedge = ptaxisinvmass->GetBinUpEdge(k);
957
958 ((TAxis *)invmassss->GetAxis(0))->SetRange(k,k);
959 ((TAxis *)invmassr->GetAxis(0))->SetRange(k,k);
960 ((TAxis *)invmassos->GetAxis(0))->SetRange(k,k);
961
962 invmassosptproj[k-1] = invmassos->Projection(2);
963 invmassssptproj[k-1] = invmassss->Projection(2);
964 invmassrptproj[k-1] = invmassr->Projection(2);
965 invmassgammaptproj[k-1] = 0x0;
966 invmasspi0ptproj[k-1] = 0x0;
967 invmassetaptproj[k-1] = 0x0;
968 invmassCptproj[k-1] = 0x0;
969 invmassBptproj[k-1] = 0x0;
970 if(invmassgamma) invmassgammaptproj[k-1] = invmassgamma->Projection(2);
971 if(invmasspi0) invmasspi0ptproj[k-1] = invmasspi0->Projection(2);
972 if(invmasseta) invmassetaptproj[k-1] = invmasseta->Projection(2);
973 if(invmassC) invmassCptproj[k-1] = invmassC->Projection(2);
974 if(invmassB) invmassBptproj[k-1] = invmassB->Projection(2);
975
976 invmassdiffptproj[k-1] = (TH1D *) invmassosptproj[k-1]->Clone();
977 TString name("Invmassdiffptbin");
978 name += k;
979 invmassdiffptproj[k-1]->SetName(name);
980 invmassdiffptproj[k-1]->Add(invmassssptproj[k-1],-1.0);
981
982 TString namee("p_{T} tagged from ");
983 namee += lowedge;
984 namee += " GeV/c to ";
985 namee += upedge;
986 namee += " GeV/c";
987
988 invmassosptproj[k-1]->SetTitle((const char*)namee);
989 invmassssptproj[k-1]->SetTitle((const char*)namee);
990 invmassrptproj[k-1]->SetTitle((const char*)namee);
991 invmassdiffptproj[k-1]->SetTitle((const char*)namee);
992 if(invmassgammaptproj[k-1]) invmassgammaptproj[k-1]->SetTitle((const char*)namee);
993 if(invmasspi0ptproj[k-1]) invmasspi0ptproj[k-1]->SetTitle((const char*)namee);
994 if(invmassetaptproj[k-1]) invmassetaptproj[k-1]->SetTitle((const char*)namee);
995 if(invmassCptproj[k-1]) invmassCptproj[k-1]->SetTitle((const char*)namee);
996 if(invmassBptproj[k-1]) invmassBptproj[k-1]->SetTitle((const char*)namee);
997
998
999
1000 invmassosptproj[k-1]->SetStats(0);
1001 invmassssptproj[k-1]->SetStats(0);
1002 invmassrptproj[k-1]->SetStats(0);
1003 invmassdiffptproj[k-1]->SetStats(0);
1004 if(invmassgammaptproj[k-1]) invmassgammaptproj[k-1]->SetStats(0);
1005 if(invmasspi0ptproj[k-1]) invmasspi0ptproj[k-1]->SetStats(0);
1006 if(invmassetaptproj[k-1]) invmassetaptproj[k-1]->SetStats(0);
1007 if(invmassCptproj[k-1]) invmassCptproj[k-1]->SetStats(0);
1008 if(invmassBptproj[k-1]) invmassBptproj[k-1]->SetStats(0);
1009
1010 Double_t yieldf = invmassdiffptproj[k-1]->Integral();
1011 if(invmassetaptproj[k-1] && invmasspi0ptproj[k-1] && invmassgammaptproj[k-1] && invmassCptproj[k-1] && invmassBptproj[k-1]) {
1012 Double_t yieldg = invmassetaptproj[k-1]->Integral() + invmasspi0ptproj[k-1]->Integral() + invmassgammaptproj[k-1]->Integral();
1013 yieldPtSourcesMC->SetBinContent(k,yieldg);
1014
1015 Double_t yieldsignal = invmassCptproj[k-1]->Integral() + invmassBptproj[k-1]->Integral();
1016 yieldPtSignalCutMC->SetBinContent(k,yieldsignal);
1017 }
1018
1019 yieldPtFound->SetBinContent(k,yieldf);
1020
1021 canvas->cd(k);
1022 invmassosptproj[k-1]->Draw();
1023 invmassssptproj[k-1]->Draw("same");
1024 invmassdiffptproj[k-1]->Draw("same");
1025 invmassrptproj[k-1]->Draw("same");
1026 TLegend *legiv = new TLegend(0.4,0.6,0.89,0.89);
1027 legiv->AddEntry(invmassosptproj[k-1],"Opposite signs","p");
1028 legiv->AddEntry(invmassssptproj[k-1],"Same signs","p");
1029 legiv->AddEntry(invmassdiffptproj[k-1],"(Opposite - Same) signs","p");
1030 legiv->AddEntry(invmassrptproj[k-1],"rotated","p");
1031 if(invmassgammaptproj[k-1]) legiv->AddEntry(invmassgammaptproj[k-1],"e^{+}e^{-} from #gamma","p");
1032 if(invmasspi0ptproj[k-1]) legiv->AddEntry(invmasspi0ptproj[k-1],"e^{+}e^{-} from #pi^{0}","p");
1033 if(invmassetaptproj[k-1]) legiv->AddEntry(invmassetaptproj[k-1],"e^{+}e^{-} from #eta","p");
1034 legiv->Draw("same");
1035
1036 }
1037
1038 yieldPtFound->SetStats(0);
1039 if(yieldPtSourcesMC) yieldPtSourcesMC->SetStats(0);
1040 if(yieldPtSignalCutMC) yieldPtSignalCutMC->SetStats(0);
1041
1042 TCanvas * canvasfin =new TCanvas("results","results",800,800);
1043 canvasfin->cd(1);
1044 yieldPtFound->Draw();
1045 if(yieldPtSourcesMC && yieldPtSignalCutMC) {
1046 yieldPtSourcesMC->Draw("same");
1047 yieldPtSignalCutMC->Draw("same");
1048 TLegend *lega = new TLegend(0.4,0.6,0.89,0.89);
1049 lega->AddEntry(yieldPtFound,"Contributions found","p");
1050 lega->AddEntry(yieldPtSourcesMC,"Contributions of e^{+}e^{-} from #gamma, #pi^{0} and #eta","p");
1051 lega->AddEntry(yieldPtSignalCutMC,"Contributions of e^{+}e^{-} from C and B","p");
1052 lega->Draw("same");
1053 }
1054
1055
1056
1057 if(!fListPostProcess) fListPostProcess = new TList();
1058 fListPostProcess->SetName("ListPostProcess");
1059
1060 for(Int_t k=0; k < nbinsptinvmass; k++){
1061 fListPostProcess->AddAt(invmassosptproj[k],kOos+kNOutput*k);
1062 fListPostProcess->AddAt(invmassssptproj[k],kOss+kNOutput*k);
1063 fListPostProcess->AddAt(invmassrptproj[k],kOr+kNOutput*k);
1064 fListPostProcess->AddAt(invmassdiffptproj[k],kOdiff+kNOutput*k);
1065 if(invmassgammaptproj[k]) fListPostProcess->AddAt(invmassgammaptproj[k],kNOutput*nbinsptinvmass+kNMCInfo*k+kElectronFromGamma);
1066 if(invmasspi0ptproj[k]) fListPostProcess->AddAt(invmasspi0ptproj[k],kNOutput*nbinsptinvmass+kNMCInfo*k+kElectronFromPi0);
1067 if(invmassetaptproj[k]) fListPostProcess->AddAt(invmassetaptproj[k],kNOutput*nbinsptinvmass+kNMCInfo*k+kElectronFromEta);
1068 if(invmassCptproj[k]) fListPostProcess->AddAt(invmassCptproj[k],kNOutput*nbinsptinvmass+kNMCInfo*k+kElectronFromC);
1069 if(invmassBptproj[k]) fListPostProcess->AddAt(invmassBptproj[k],kNOutput*nbinsptinvmass+kNMCInfo*k+kElectronFromB);
1070 }
1071
1072 fListPostProcess->AddAt(yieldPtFound,kNOutput*nbinsptinvmass+kNMCInfo*nbinsptinvmass);
1073 if(yieldPtSourcesMC) fListPostProcess->AddAt(yieldPtSourcesMC,kNOutput*nbinsptinvmass+kNMCInfo*nbinsptinvmass+1);
1074 if(yieldPtSignalCutMC) fListPostProcess->AddAt(yieldPtSignalCutMC,kNOutput*nbinsptinvmass+kNMCInfo*nbinsptinvmass+2);
1075
1076 // delete dynamic array
1077 delete[] invmassosptproj;
1078 delete[] invmassssptproj;
1079 delete[] invmassrptproj;
1080 delete[] invmassdiffptproj;
1081 delete[] invmassgammaptproj;
1082 delete[] invmasspi0ptproj;
1083 delete[] invmassetaptproj;
1084 delete[] invmassCptproj;
1085 delete[] invmassBptproj;
1086
1087}
1088//_______________________________________________________________________________________________
1089void AliHFEelecbackground::Plot() const
1090{
1091 //
1092 // Plot the output
1093 //
1094
1095 if(!fList) return;
1096
1097 // opening angle
1098 THnSparseF *openinganglepp = dynamic_cast<THnSparseF *>(fList->FindObject("openinganglepp"));
1099 THnSparseF *openinganglenn = dynamic_cast<THnSparseF *>(fList->FindObject("openinganglenn"));
1100 THnSparseF *openingangless = dynamic_cast<THnSparseF *>(fList->FindObject("openingangless"));
1101 THnSparseF *openingangler = dynamic_cast<THnSparseF *>(fList->FindObject("openingangler"));
1102 THnSparseF *openingangleos = dynamic_cast<THnSparseF *>(fList->FindObject("openingangleos"));
1103
1104 THnSparseF *openinganglegamma = dynamic_cast<THnSparseF *>(fList->FindObject("openinganglegamma"));
1105 THnSparseF *openinganglepi0 = dynamic_cast<THnSparseF *>(fList->FindObject("openinganglepi0"));
1106 THnSparseF *openingangleC = dynamic_cast<THnSparseF *>(fList->FindObject("openingangleC"));
1107 THnSparseF *openingangleB = dynamic_cast<THnSparseF *>(fList->FindObject("openingangleB"));
1108 THnSparseF *openingangleeta = dynamic_cast<THnSparseF *>(fList->FindObject("openingangleeta"));
1109
1110
1111 // invariant mass
1112 THnSparseF *invmasspp = dynamic_cast<THnSparseF *>(fList->FindObject("invmasspp"));
1113 THnSparseF *invmassnn = dynamic_cast<THnSparseF *>(fList->FindObject("invmassnn"));
1114 THnSparseF *invmassss = dynamic_cast<THnSparseF *>(fList->FindObject("invmassss"));
1115 THnSparseF *invmassr = dynamic_cast<THnSparseF *>(fList->FindObject("invmassr"));
1116 THnSparseF *invmassos = dynamic_cast<THnSparseF *>(fList->FindObject("invmassos"));
1117
1118 THnSparseF *invmassgamma = dynamic_cast<THnSparseF *>(fList->FindObject("invmassgamma"));
1119 THnSparseF *invmasspi0 = dynamic_cast<THnSparseF *>(fList->FindObject("invmasspi0"));
1120 THnSparseF *invmassC = dynamic_cast<THnSparseF *>(fList->FindObject("invmassC"));
1121 THnSparseF *invmassB = dynamic_cast<THnSparseF *>(fList->FindObject("invmassB"));
1122 THnSparseF *invmasseta = dynamic_cast<THnSparseF *>(fList->FindObject("invmasseta"));
1123
1124 // Projection over all pt
1125 TH1D *openingangleppproj = openinganglepp->Projection(2);
1126 TH1D *openinganglennproj = openinganglenn->Projection(2);
1127 TH1D *openinganglessproj = openingangless->Projection(2);
1128 TH1D *openinganglerproj = openingangler->Projection(2);
1129 TH1D *openingangleosproj = openingangleos->Projection(2);
1130
1131 TH1D *openinganglegammaproj = 0x0;
1132 TH1D *openinganglepi0proj = 0x0;
1133 TH1D *openingangleCproj = 0x0;
1134 TH1D *openingangleBproj = 0x0;
1135 TH1D *openingangleetaproj = 0x0;
1136 if(openinganglegamma) openinganglegammaproj = openinganglegamma->Projection(2);
1137 if(openinganglepi0) openinganglepi0proj = openinganglepi0->Projection(2);
1138 if(openingangleC) openingangleCproj = openingangleC->Projection(2);
1139 if(openingangleB) openingangleBproj = openingangleB->Projection(2);
1140 if(openingangleeta) openingangleetaproj = openingangleeta->Projection(2);
1141
1142
1143 TH1D *invmassppproj = invmasspp->Projection(2);
1144 TH1D *invmassnnproj = invmassnn->Projection(2);
1145 TH1D *invmassssproj = invmassss->Projection(2);
1146 TH1D *invmassrproj = invmassr->Projection(2);
1147 TH1D *invmassosproj = invmassos->Projection(2);
1148
1149 TH1D *invmassgammaproj = 0x0;
1150 TH1D *invmasspi0proj = 0x0;
1151 TH1D *invmassCproj = 0x0;
1152 TH1D *invmassBproj = 0x0;
1153 TH1D *invmassetaproj = 0x0;
1154 if(invmassgamma) invmassgammaproj = invmassgamma->Projection(2);
1155 if(invmasspi0) invmasspi0proj = invmasspi0->Projection(2);
1156 if(invmassC) invmassCproj = invmassC->Projection(2);
1157 if(invmassB) invmassBproj = invmassB->Projection(2);
1158 if(invmasseta) invmassetaproj = invmasseta->Projection(2);
1159
1160 openingangleppproj->SetStats(0);
1161 openinganglennproj->SetStats(0);
1162 openinganglessproj->SetStats(0);
1163 openinganglerproj->SetStats(0);
1164 openingangleosproj->SetStats(0);
1165 if(openinganglegammaproj) openinganglegammaproj->SetStats(0);
1166 if(openinganglepi0proj) openinganglepi0proj->SetStats(0);
1167 if(openingangleCproj) openingangleCproj->SetStats(0);
1168 if(openingangleBproj) openingangleBproj->SetStats(0);
1169 if(openingangleetaproj) openingangleetaproj->SetStats(0);
1170
1171 invmassppproj->SetStats(0);
1172 invmassnnproj->SetStats(0);
1173 invmassssproj->SetStats(0);
1174 invmassrproj->SetStats(0);
1175 invmassosproj->SetStats(0);
1176 if(invmassgammaproj) invmassgammaproj->SetStats(0);
1177 if(invmasspi0proj) invmasspi0proj->SetStats(0);
1178 if(invmassCproj) invmassCproj->SetStats(0);
1179 if(invmassBproj) invmassBproj->SetStats(0);
1180 if(invmassetaproj) invmassetaproj->SetStats(0);
1181
1182 openingangleppproj->SetTitle("");
1183 openinganglennproj->SetTitle("");
1184 openinganglessproj->SetTitle("");
1185 openinganglerproj->SetTitle("");
1186 openingangleosproj->SetTitle("");
1187 if(openinganglegammaproj) openinganglegammaproj->SetTitle("");
1188 if(openinganglepi0proj) openinganglepi0proj->SetTitle("");
1189 if(openingangleCproj) openingangleCproj->SetTitle("");
1190 if(openingangleBproj) openingangleBproj->SetTitle("");
1191 if(openingangleetaproj) openingangleetaproj->SetTitle("");
1192
1193 invmassppproj->SetTitle("");
1194 invmassnnproj->SetTitle("");
1195 invmassssproj->SetTitle("");
1196 invmassrproj->SetTitle("");
1197 invmassosproj->SetTitle("");
1198 if(invmassgammaproj) invmassgammaproj->SetTitle("");
1199 if(invmasspi0proj) invmasspi0proj->SetTitle("");
1200 if(invmassCproj) invmassCproj->SetTitle("");
1201 if(invmassBproj) invmassBproj->SetTitle("");
1202 if(invmassetaproj) invmassetaproj->SetTitle("");
1203
1204 // Draw histograms for opening angle
1205 TCanvas * copeningangle =new TCanvas("openingangle","Openingangle",800,800);
1206 copeningangle->cd();
1207 openingangleppproj->Draw();
1208 openinganglennproj->Draw("same");
1209 openinganglessproj->Draw("same");
1210 openinganglerproj->Draw("same");
1211 openingangleosproj->Draw("same");
1212 if(openinganglegammaproj) openinganglegammaproj->Draw("same");
1213 if(openinganglepi0proj) openinganglepi0proj->Draw("same");
1214 if(openingangleCproj) openingangleCproj->Draw("same");
1215 if(openingangleBproj) openingangleBproj->Draw("same");
1216 if(openingangleetaproj) openingangleetaproj->Draw("same");
1217 TLegend *lego = new TLegend(0.4,0.6,0.89,0.89);
1218 lego->AddEntry(openingangleppproj,"positive-positive","p");
1219 lego->AddEntry(openinganglennproj,"negative-negative","p");
1220 lego->AddEntry(openinganglessproj,"same-sign","p");
1221 lego->AddEntry(openinganglerproj,"rotated","p");
1222 lego->AddEntry(openingangleosproj,"positive-negative","p");
1223 if(openinganglegammaproj) lego->AddEntry(openinganglegammaproj,"e^{+}e^{-} from #gamma","p");
1224 if(openinganglepi0proj) lego->AddEntry(openinganglepi0proj,"e^{+}e^{-} from #pi^{0}","p");
1225 if(openingangleCproj) lego->AddEntry(openingangleCproj,"e^{+}e^{-} from c","p");
1226 if(openingangleBproj) lego->AddEntry(openingangleBproj,"e^{+}e^{-} from b","p");
1227 if(openingangleetaproj) lego->AddEntry(openingangleetaproj,"e^{+}e^{-} from #eta","p");
1228 lego->Draw("same");
1229
1230 // Draw histograms for opening angle
1231 TCanvas * cinvmass =new TCanvas("invmass","Invmass",800,800);
1232 cinvmass->cd();
1233 invmassppproj->Draw();
1234 invmassnnproj->Draw("same");
1235 invmassssproj->Draw("same");
1236 invmassrproj->Draw("same");
1237 invmassosproj->Draw("same");
1238 if(invmassgammaproj) invmassgammaproj->Draw("same");
1239 if(invmasspi0proj) invmasspi0proj->Draw("same");
1240 if(invmassCproj) invmassCproj->Draw("same");
1241 if(invmassBproj) invmassBproj->Draw("same");
1242 if(invmassetaproj) invmassetaproj->Draw("same");
1243 TLegend *legi = new TLegend(0.4,0.6,0.89,0.89);
1244 legi->AddEntry(invmassppproj,"positive-positive","p");
1245 legi->AddEntry(invmassnnproj,"negative-negative","p");
1246 legi->AddEntry(invmassssproj,"same-sign","p");
1247 legi->AddEntry(invmassrproj,"rotated","p");
1248 legi->AddEntry(invmassosproj,"positive-negative","p");
1249 if(invmassgammaproj) legi->AddEntry(invmassgammaproj,"e^{+}e^{-} from #gamma","p");
1250 if(invmasspi0proj) legi->AddEntry(invmasspi0proj,"e^{+}e^{-} from #pi^{0}","p");
1251 if(invmassCproj) legi->AddEntry(invmassCproj,"e^{+}e^{-} from c","p");
1252 if(invmassBproj) legi->AddEntry(invmassBproj,"e^{+}e^{-} from b","p");
1253 if(invmassetaproj) legi->AddEntry(invmassetaproj,"e^{+}e^{-} from #eta","p");
1254 legi->Draw("same");
1255
1256}