]> 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
Updates, adding EMCal (committing on behalf of Silvia)
[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 **************************************************************************/
27de2dfb 15
16/* $Id$ */
17
02524e30 18//
19//
20// First implementation of a class
21// to reject tagged electron coming from conversion, pi0 and eta
22// by calculating the e+e- invariant mass
23// of the tagged electron with other tracks
24// after looser cuts for the partner.
25// PostProcess should extract the background yield
26// If running with MC, it can be compared to the expected background yield
27//
28// Authors:
29// Raphaelle Bailhache <rbailhache@ikf.uni-frankfurt.de > <R.Bailhache@gsi.de >
30//
31//
32
33#include <THnSparse.h>
34#include <TParticle.h>
35#include <TFile.h>
36#include <TCanvas.h>
37#include <TAxis.h>
38#include <TH1D.h>
d2af20c5 39#include <TH2D.h>
70da6c5a 40#include <TH1F.h>
41#include <TH2F.h>
42#include <TF1.h>
02524e30 43#include <TString.h>
44#include <TLegend.h>
70da6c5a 45#include <TStyle.h>
02524e30 46
47#include <AliESDEvent.h>
48#include <AliAODEvent.h>
49#include <AliESDtrack.h>
50#include <AliAODTrack.h>
51#include "AliHFEelecbackground.h"
52#include <AliMCEvent.h>
53#include <AliStack.h>
70da6c5a 54#include <AliKFParticle.h>
55#include "AliCFContainer.h"
56#include "AliHFEpid.h"
57#include "AliESDpid.h"
69ac0e6f 58#include "AliLog.h"
70da6c5a 59#include "AliITSPIDResponse.h"
60#include "AliTPCPIDResponse.h"
02524e30 61
62ClassImp(AliHFEelecbackground)
63
d2af20c5 64Bool_t AliHFEelecbackground::fgUseMCPID = kFALSE;
65const Double_t AliHFEelecbackground::fgkMe = 0.00051099892;
02524e30 66
67//___________________________________________________________________________________________
68AliHFEelecbackground::AliHFEelecbackground():
69ac0e6f 69 fhtmp(0x0)
70 ,fhtmpf(0x0)
71 ,fhtmpp(0x0)
72 ,fESD1(0x0)
02524e30 73 ,fAOD1(0x0)
74 ,fMCEvent(0x0)
75 ,fBz(0)
76 ,fkVertex(0x0)
02524e30 77 ,fPtESD(0.0)
78 ,fIndexTrack(0)
02524e30 79 ,fPdg(0)
80 ,fLabMother(-1)
81 ,fIsFrom(-1)
82 ,fMotherGamma(-1)
83 ,fMotherPi0(-1)
84 ,fMotherC(-1)
85 ,fMotherB(-1)
86 ,fMotherEta(-1)
87 ,fIsPartner(kFALSE)
d2af20c5 88 ,fIsSplittedTrack(kFALSE)
70da6c5a 89 ,fOpeningAngleCut(0.35)
90 ,fInvMassCut(140.0)
91 ,fChi2NdfCut(999999999.0)
92 ,fUseAliKFCode(kTRUE)
93 ,fSharedClusterCut(kFALSE)
94 ,fRequireITSStandalone(0)
95 ,fMinNbCls(2)
96 ,fMinITSChi2(10.0)
97 ,fMinNbClsSDDSPD(2)
98 ,fPIDPartner(kFALSE)
99 ,fPIDMethodPartner(0x0)
100 ,fPIDMethodPartnerITS(0x0)
101 ,fDebugLevel(0)
02524e30 102 ,fList(0x0)
103 ,fListPostProcess(0x0)
104{
105 //
106 // Default constructor
107 //
69ac0e6f 108 for(Int_t k =0; k < 10; k++) {
109 fCuts[k] = kFALSE;
110 }
02524e30 111
112}
113
114//_______________________________________________________________________________________________
115AliHFEelecbackground::AliHFEelecbackground(const AliHFEelecbackground &p):
116 TObject(p)
69ac0e6f 117 ,fhtmp(0x0)
118 ,fhtmpf(0x0)
119 ,fhtmpp(0x0)
02524e30 120 ,fESD1(0x0)
121 ,fAOD1(0x0)
122 ,fMCEvent(0x0)
123 ,fBz(p.fBz)
124 ,fkVertex(p.fkVertex)
02524e30 125 ,fPtESD(p.fPtESD)
126 ,fIndexTrack(0)
02524e30 127 ,fPdg(0)
128 ,fLabMother(-1)
129 ,fIsFrom(-1)
130 ,fMotherGamma(-1)
131 ,fMotherPi0(-1)
132 ,fMotherC(-1)
133 ,fMotherB(-1)
134 ,fMotherEta(-1)
135 ,fIsPartner(kFALSE)
d2af20c5 136 ,fIsSplittedTrack(kFALSE)
70da6c5a 137 ,fOpeningAngleCut(0.35)
138 ,fInvMassCut(140.0)
139 ,fChi2NdfCut(999999999.0)
140 ,fUseAliKFCode(kTRUE)
141 ,fSharedClusterCut(kFALSE)
142 ,fRequireITSStandalone(0)
143 ,fMinNbCls(2)
144 ,fMinITSChi2(10.0)
145 ,fMinNbClsSDDSPD(2)
146 ,fPIDPartner(kFALSE)
147 ,fPIDMethodPartner(0x0)
148 ,fPIDMethodPartnerITS(0x0)
149 ,fDebugLevel(0)
02524e30 150 ,fList(0x0)
151 ,fListPostProcess(0x0)
152{
153 //
154 // Copy constructor
155 //
69ac0e6f 156 for(Int_t k =0; k < 10; k++) {
157 fCuts[k] = kFALSE;
158 }
02524e30 159}
160
161//_______________________________________________________________________________________________
162AliHFEelecbackground&
163AliHFEelecbackground::operator=(const AliHFEelecbackground &)
164{
165 //
166 // Assignment operator
167 //
168
169 AliInfo("Not yet implemented.");
170 return *this;
171}
172
173//_______________________________________________________________________________________________
174AliHFEelecbackground::~AliHFEelecbackground()
175{
176 //
177 // Destructor
178 //
70da6c5a 179 if(fPIDMethodPartner) delete fPIDMethodPartner;
180 if(fPIDMethodPartnerITS) delete fPIDMethodPartnerITS;
02524e30 181
02524e30 182 if(fListPostProcess){
3a72645a 183 fListPostProcess->SetOwner(kTRUE);
02524e30 184 delete fListPostProcess;
185 }
3a72645a 186
187/*
188 if(fhtmp) delete fhtmp;
189 if(fhtmpf) delete fhtmpf;
190 if(fhtmpp) delete fhtmpp;
191*/
192
02524e30 193}
194//___________________________________________________________________________________________
70da6c5a 195Bool_t AliHFEelecbackground::Load(const Char_t * filename)
02524e30 196{
197 //
198 // Generic container loader
199 //
200
201 if(!TFile::Open(filename)){
202 return kFALSE;
203 }
204 TList *o = 0x0;
205 if(!(o = (TList*)gFile->Get("Results"))){
206 return kFALSE;
207 }
208 TList *oe = 0x0;
209 if(!(oe = (TList*)dynamic_cast<TList *>(o->FindObject("HFEelecbackground")))){
210 return kFALSE;
211 }
212 fList = (TList*)oe->Clone("HFEelecbackground");
213 gFile->Close();
214 return kTRUE;
215}
70da6c5a 216//___________________________________________________________________________________________
217Bool_t AliHFEelecbackground::Load(TList * const outputlist)
218{
219 //
220 // Generic container loader
221 //
222 if(!outputlist) return kFALSE;
223 else fList = (TList*)outputlist->Clone("HFEelecbackground");
224 return kTRUE;
225}
02524e30 226//_______________________________________________________________________________________________
d2af20c5 227void AliHFEelecbackground::Reset()
228{
229 //
230 // Reset variables
231 //
232 fPtESD = 0.0;
233 fIndexTrack = -1;
234 fPdg = -1;
235 fLabMother = -1;
236 fIsFrom = -1;
237 fMotherGamma = -1;
238 fMotherPi0 = -1;
239 fMotherC = -1;
240 fMotherB = -1;
241 fMotherEta = -1;
242 fIsPartner = kFALSE;
243 fIsSplittedTrack = kFALSE;
70da6c5a 244 for(Int_t id = 0; id < 10; id++) {
245 fCuts[id] = kFALSE;
246 }
d2af20c5 247
248}
249//_______________________________________________________________________________________________
70da6c5a 250void AliHFEelecbackground::CreateHistograms(TList * const qaList)
02524e30 251{
252 //
253 // create histograms
254 //
255 if(!qaList) return;
70da6c5a 256
02524e30 257 fList = qaList;
70da6c5a 258 fList->SetName("HFEelecbackground");
259
260 //////////
02524e30 261 // bins
70da6c5a 262 /////////
02524e30 263
67fe7bd0 264 const Int_t nBinsPt = 25;
70da6c5a 265 Double_t minPt = 0.01;
02524e30 266 Double_t maxPt = 10.0;
267
67fe7bd0 268 const Int_t nBinsPtMore = 100;
70da6c5a 269 Double_t minPtMore = 0.01;
270 Double_t maxPtMore = 10.0;
271
67fe7bd0 272 const Int_t nBinsInv = 50;
02524e30 273 Double_t minInv = 0.0;
274 Double_t maxInv = 0.2;
275
67fe7bd0 276 const Int_t nBinsOp = 50;
02524e30 277 Double_t minOp = 0.0;
278 Double_t maxOp = 2;
279
67fe7bd0 280 const Int_t nBinsCh = 4;
70da6c5a 281 Double_t minCh = 0.0;
282 Double_t maxCh = 4.0;
283
67fe7bd0 284 Double_t binLimLogPt[nBinsPt+1];
285 Double_t binLimPt[nBinsPt+1];
02524e30 286 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 ;
287 for(Int_t i=0; i<=nBinsPt; i++) binLimPt[i]=(Double_t)TMath::Power(10,binLimLogPt[i]);
288
67fe7bd0 289 Double_t binLimLogPtMore[nBinsPtMore+1];
290 Double_t binLimPtMore[nBinsPtMore+1];
70da6c5a 291 for(Int_t i=0; i<=nBinsPtMore; i++) binLimLogPtMore[i]=(Double_t)TMath::Log10(minPtMore) + (TMath::Log10(maxPtMore)-TMath::Log10(minPtMore))/nBinsPtMore*(Double_t)i ;
292 for(Int_t i=0; i<=nBinsPtMore; i++) binLimPtMore[i]=(Double_t)TMath::Power(10,binLimLogPtMore[i]);
293
67fe7bd0 294 Double_t binLimInv[nBinsInv+1];
02524e30 295 for(Int_t i=0; i<=nBinsInv; i++) binLimInv[i]=(Double_t)minInv + (maxInv-minInv) /nBinsInv*(Double_t)i ;
296
67fe7bd0 297 Double_t binLimOp[nBinsOp+1];
02524e30 298 for(Int_t i=0; i<=nBinsOp; i++) binLimOp[i]=(Double_t)minOp + (maxOp-minOp) /nBinsOp*(Double_t)i ;
299
67fe7bd0 300 Double_t binLimCh[nBinsCh+1];
70da6c5a 301 for(Int_t i=0; i<=nBinsCh; i++) binLimCh[i]=(Double_t)minCh + (maxCh-minCh) /nBinsCh*(Double_t)i ;
02524e30 302
70da6c5a 303 const Int_t nvarData = 5;
304 // pt reconstructed tagged e
305 // pt reconstructed mother
306 // opening angle
307 // invariant mass
308 // Data: charge (0-opposite sign, 1-like sign ++, 2-like sign --, 3-rotated tracks)
02524e30 309
70da6c5a 310 Int_t iBinData[nvarData];
311 iBinData[0]=nBinsPt;
312 iBinData[1]=nBinsPt;
313 iBinData[2]=nBinsOp;
314 iBinData[3]=nBinsInv;
315 iBinData[4]=nBinsCh;
02524e30 316
02524e30 317 //
70da6c5a 318 // Opening angle and invariant mass
02524e30 319 //
320
70da6c5a 321 THnSparseF *hsSparseData = new THnSparseF("OpeningangleinvmassData","",nvarData,iBinData);
67fe7bd0 322 hsSparseData->SetBinEdges(0,&binLimPt[0]);
323 hsSparseData->SetBinEdges(1,&binLimPt[0]);
324 hsSparseData->SetBinEdges(2,&binLimOp[0]);
325 hsSparseData->SetBinEdges(3,&binLimInv[0]);
326 hsSparseData->SetBinEdges(4,&binLimCh[0]);
70da6c5a 327 hsSparseData->Sumw2();
02524e30 328
70da6c5a 329 fList->AddAt(hsSparseData,kDatai);
02524e30 330
70da6c5a 331 //
332 // Radius, DCA and Chi2Ndf
333 //
02524e30 334
70da6c5a 335 TH1F *dataRadiusHisto = new TH1F("DataRadius","", 200, 0.0, 200.0); // recontructed radius from the AliKF of the e+e- pair
336 fList->AddAt(dataRadiusHisto,kDatar);
d2af20c5 337
70da6c5a 338 TH1F *dataDcaHisto = new TH1F("DataDCA","", 100, 0.0, 6.0); // dca distribution
339 fList->AddAt(dataDcaHisto,kDatadca);
02524e30 340
70da6c5a 341 TH1F *dataChi2NdfHisto = new TH1F("DataChi2Ndf","", 100, 0.0, 5.0); // chi2Ndf distribution
342 fList->AddAt(dataChi2NdfHisto,kDatachi2Ndf);
02524e30 343
70da6c5a 344
02524e30 345 if(HasMCData()) {
70da6c5a 346
347 //
348 // Opening angle and invariant mass with MC infos
349 //
350
351 const Int_t nvarMCo = 6;
352 // pt reconstructed tagged e
353 // pt reconstructed mother
354 // opening angle
355 // invariant mass
356 // MC: 0-FromBackground, 1-FromGamma, 2-FromPi0, 3-FromEta, 4-FromC, 5-FromB
357 // MCSplitted: 0-not, 1-splittedOs, 2-ksplittedSs
02524e30 358
02524e30 359
67fe7bd0 360 const Int_t nBinsMCOrigin = 6;
70da6c5a 361 Double_t minMCOrigin = 0.0;
362 Double_t maxMCOrigin = 6.0;
363
67fe7bd0 364 Double_t binLimMCOrigin[nBinsMCOrigin+1];
70da6c5a 365 for(Int_t i=0; i<=nBinsMCOrigin; i++) binLimMCOrigin[i]=(Double_t)minMCOrigin + (maxMCOrigin-minMCOrigin) /nBinsMCOrigin*(Double_t)i ;
d2af20c5 366
67fe7bd0 367 const Int_t nBinsMCSplitted = 3;
70da6c5a 368 Double_t minMCSplitted = 0.0;
369 Double_t maxMCSplitted = 3.0;
370
67fe7bd0 371 Double_t binLimMCSplitted[nBinsMCSplitted+1];
70da6c5a 372 for(Int_t i=0; i<=nBinsMCSplitted; i++) binLimMCSplitted[i]=(Double_t)minMCSplitted + (maxMCSplitted-minMCSplitted) /nBinsMCSplitted*(Double_t)i ;
373
374 Int_t iBinMCo[nvarMCo];
375 iBinMCo[0]=nBinsPt;
376 iBinMCo[1]=nBinsPt;
377 iBinMCo[2]=nBinsOp;
378 iBinMCo[3]=nBinsInv;
379 iBinMCo[4]=nBinsMCOrigin;
380 iBinMCo[5]=nBinsMCSplitted;
381
382 THnSparseF *hsSparseMCo = new THnSparseF("OpeningangleinvmassMC","",nvarMCo,iBinMCo);
67fe7bd0 383 hsSparseMCo->SetBinEdges(0,&binLimPt[0]);
384 hsSparseMCo->SetBinEdges(1,&binLimPt[0]);
385 hsSparseMCo->SetBinEdges(2,&binLimOp[0]);
386 hsSparseMCo->SetBinEdges(3,&binLimInv[0]);
387 hsSparseMCo->SetBinEdges(4,&binLimMCOrigin[0]);
388 hsSparseMCo->SetBinEdges(5,&binLimMCSplitted[0]);
70da6c5a 389 hsSparseMCo->Sumw2();
390
391 fList->AddAt(hsSparseMCo,kMCo);
392
393 //
394 // Radius, DCA and Chi2Ndf with MC info
395 //
396
397 TH2F *mcRadiusHisto = new TH2F("MCRadius","", 200, 0.0, 200.0,6,-0.5,5.5); // recontructed radius from the AliKF of the e+e- pair
398 fList->AddAt(mcRadiusHisto,kMCr);
399
400 TH2F *mcDcaHisto = new TH2F("MCDCA","", 100, 0.0, 6.0,6,-0.5,5.5); // dca distribution
401 fList->AddAt(mcDcaHisto,kMCdca);
402
403 TH2F *mcChi2NdfHisto = new TH2F("MCChi2Ndf","", 100, 0.0, 5.0,6,-0.5,5.5); // chi2Ndf distribution
404 fList->AddAt(mcChi2NdfHisto,kMCchi2Ndf);
405
406 //////////////////////////////////////////////////////////
407 // if fDebugLevel 1: Rejection efficiencies of the cuts
408 //////////////////////////////////////////////////////////
409
410 if(fDebugLevel > 0) {
411
412 if(HasMCData()) {
413
414 const Int_t nvarMCe = 3;
415 // pt reconstructed tagged e
416 // cut passed: 0-all, 1-Partner tracked, 2-Opposite-sign, 3-SingleTrackCutPart, 4-ShareCluster, 5-PID, 6-DCA, 7-chi2Ndf AliKF, 8-Openingangle, 9-Invmass
417 // MC: 0-FromBackground, 1-FromGamma, 2-FromPi0, 3-FromEta, 4-FromC, 5-FromB
418
67fe7bd0 419 const Int_t nBinsMCCutPassed = 10;
70da6c5a 420 Double_t minMCCutPassed = -0.5;
421 Double_t maxMCCutPassed = 9.5;
422
67fe7bd0 423 Double_t binLimMCCutPassed[nBinsMCCutPassed+1];
70da6c5a 424 for(Int_t i=0; i<=nBinsMCCutPassed; i++) binLimMCCutPassed[i]=(Double_t)minMCCutPassed + (maxMCCutPassed-minMCCutPassed) /nBinsMCCutPassed*(Double_t)i ;
425
426 Int_t iBinMCe[nvarMCe];
427 iBinMCe[0]=nBinsPt;
428 iBinMCe[1]=nBinsMCCutPassed;
429 iBinMCe[2]=nBinsMCOrigin;
430
431 THnSparseF *hsSparseMCe = new THnSparseF("CutPassedMC","",nvarMCe,iBinMCe);
67fe7bd0 432 hsSparseMCe->SetBinEdges(0,&binLimPt[0]);
433 hsSparseMCe->SetBinEdges(1,&binLimMCCutPassed[0]);
434 hsSparseMCe->SetBinEdges(2,&binLimMCOrigin[0]);
70da6c5a 435 hsSparseMCe->Sumw2();
436
437 fList->AddAt(hsSparseMCe,kMCe);
438
439 }
440 }
441
442 /////////////////////////////////////////////////////////////////
443 // if fDebugLevel 1: PIDPartCut and ShareClusters
444 /////////////////////////////////////////////////////////////////
445
446 if(fDebugLevel > 1) {
447
448 if(!fRequireITSStandalone){
449
450 //
451 // TPC
452 //
453
454 TH2F *tpcPartner0 = new TH2F("TPCPartner0","", nBinsPtMore, binLimPtMore, 200, 0.0, 700.0);
455 fList->AddAt(tpcPartner0,kMCcutPart0);
456 TH2F *tpcPartner1 = new TH2F("TPCPartner1","", nBinsPtMore, binLimPtMore, 200, 0.0, 700.0);
457 fList->AddAt(tpcPartner1,kMCcutPart1);
458
459 }
460 else {
461
462 //
463 // ITS
464 //
465
466 TH2F *itsPartner0 = new TH2F("ITSPartner0","", nBinsPtMore, binLimPtMore, 200, 0.0, 700.0);
467 fList->AddAt(itsPartner0,kMCcutPart0);
468 TH2F *itsPartner1 = new TH2F("ITSPartner1","", nBinsPtMore, binLimPtMore, 200, 0.0, 700.0);
469 fList->AddAt(itsPartner1,kMCcutPart1);
470
471 /////////////////////////////////////////////////////
472 // dEdx of the four layers for the track partner
473 /////////////////////////////////////////////////////
474 const Int_t nvarITSsignal = 5;
475
67fe7bd0 476 const Int_t nBinsITSsignal = 100;
70da6c5a 477 Double_t minITSsignal = 0.0;
478 Double_t maxITSsignal = 350.0;
479
67fe7bd0 480 Double_t binLimITSsignal[nBinsITSsignal+1];
70da6c5a 481 for(Int_t i=0; i<=nBinsITSsignal; i++) binLimITSsignal[i]=(Double_t)minITSsignal + (maxITSsignal-minITSsignal) /nBinsITSsignal*(Double_t)i ;
482
483 Int_t iBinITSsignal[nvarITSsignal];
484 iBinITSsignal[0]=nBinsPt;
485 iBinITSsignal[1]=nBinsITSsignal;
486 iBinITSsignal[2]=nBinsITSsignal;
487 iBinITSsignal[3]=nBinsITSsignal;
488 iBinITSsignal[4]=nBinsITSsignal;
489
490 THnSparseF *hsSparseITSpid = new THnSparseF("SparseITSsignal","",nvarITSsignal,iBinITSsignal);
67fe7bd0 491 hsSparseITSpid->SetBinEdges(0,&binLimPt[0]);
492 hsSparseITSpid->SetBinEdges(1,&binLimITSsignal[0]);
493 hsSparseITSpid->SetBinEdges(2,&binLimITSsignal[0]);
494 hsSparseITSpid->SetBinEdges(3,&binLimITSsignal[0]);
495 hsSparseITSpid->SetBinEdges(4,&binLimITSsignal[0]);
70da6c5a 496 hsSparseITSpid->Sumw2();
497
498 fList->AddAt(hsSparseITSpid,kMCcutPart2);
499
500 ///////////////////////////////////////////////////////////////////////////////////////
501 // dEdx of the four layers for the track partner and track to reject splitted track
502 ///////////////////////////////////////////////////////////////////////////////////////
503 const Int_t nvarITSsignalSplit = 5;
504
67fe7bd0 505 const Int_t nBinsITSSplit = 2;
70da6c5a 506 Double_t minITSSplit = 0.0;
507 Double_t maxITSSplit = 2.0;
508
67fe7bd0 509 Double_t binLimITSSplit[nBinsITSSplit+1];
70da6c5a 510 for(Int_t i=0; i<=nBinsITSSplit; i++) binLimITSSplit[i]=(Double_t)minITSSplit + (maxITSSplit-minITSSplit) /nBinsITSSplit*(Double_t)i ;
511
512
67fe7bd0 513 const Int_t nBinsITSsignalSplit = 50;
70da6c5a 514 Double_t minITSsignalSplit = -25.0;
515 Double_t maxITSsignalSplit = 25.0;
516
67fe7bd0 517 Double_t binLimITSsignalSplit[nBinsITSsignalSplit+1];
70da6c5a 518 for(Int_t i=0; i<=nBinsITSsignalSplit; i++) binLimITSsignalSplit[i]=(Double_t)minITSsignalSplit + (maxITSsignalSplit-minITSsignalSplit) /nBinsITSsignalSplit*(Double_t)i ;
519
520 Int_t iBinITSsignalSplit[nvarITSsignalSplit];
521 iBinITSsignalSplit[0]=nBinsITSSplit;
522 for(Int_t k = 1; k < 5; k++){
523 iBinITSsignalSplit[k]=nBinsITSsignalSplit;
524 }
525
526 THnSparseF *hsSparseITSpidSplit = new THnSparseF("SparseITSsignalSplit","",nvarITSsignalSplit,iBinITSsignalSplit);
67fe7bd0 527 hsSparseITSpidSplit->SetBinEdges(0,&binLimITSSplit[0]);
70da6c5a 528 for(Int_t k = 1; k < 5; k++) {
67fe7bd0 529 hsSparseITSpidSplit->SetBinEdges(k,&binLimITSsignalSplit[0]);
70da6c5a 530 }
531 hsSparseITSpidSplit->Sumw2();
532
533 fList->AddAt(hsSparseITSpidSplit,kMCcutPart3);
534
535 }
536
537 }
d2af20c5 538
02524e30 539 }
540
70da6c5a 541 //qaList->Add(fList);
542
02524e30 543}
544//_______________________________________________________________________________________________
545void AliHFEelecbackground::PairAnalysis(AliESDtrack* const track, AliESDtrack* const trackPart)
546{
547 //
548 // calculate (tagged e-partner) dca, opening angle, invariant mass
549 //
02524e30 550
70da6c5a 551 /////////////////////
552 // pt tagged
553 //////////////////////
554 TVector3 v3Dtagged;
67fe7bd0 555 Double_t pxyz[3];
556 track->PxPyPz(&pxyz[0]);
70da6c5a 557 v3Dtagged.SetXYZ(pxyz[0],pxyz[1],pxyz[2]);
558 fPtESD = TMath::Sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]);
559
02524e30 560 ////////////////////////
561 // Take label
562 ////////////////////////
563 Int_t indexTrack = TMath::Abs(track->GetLabel());
564 Int_t indexTrackPart = TMath::Abs(trackPart->GetLabel());
70da6c5a 565
02524e30 566 /////////////////////////
567 // If MC data
568 ////////////////////////
569
570 if(HasMCData()) {
571
572 // Take info track if not already done
573 if(indexTrack!= fIndexTrack) {
70da6c5a 574
575 for(Int_t id = 0; id < 10; id++) {
576 fCuts[id] = kFALSE;
577 }
d2af20c5 578
70da6c5a 579 fIsFrom = kElectronFromBackground;
d2af20c5 580
02524e30 581 fPdg = GetPdg(indexTrack);
582 fLabMother = GetLabMother(indexTrack);
583
584 fMotherGamma = IsMotherGamma(indexTrack);
585 if((fMotherGamma != -1) && ((TMath::Abs(fPdg)) == 11)) fIsFrom = kElectronFromGamma;
586 fMotherPi0 = IsMotherPi0(indexTrack);
587 if((fMotherPi0 != -1) && ((TMath::Abs(fPdg)) == 11)) fIsFrom = kElectronFromPi0;
588 fMotherC = IsMotherC(indexTrack);
589 if((fMotherC != -1) && ((TMath::Abs(fPdg)) == 11)) fIsFrom = kElectronFromC;
590 fMotherB = IsMotherB(indexTrack);
591 if((fMotherB != -1) && ((TMath::Abs(fPdg)) == 11)) fIsFrom = kElectronFromB;
592 fMotherEta = IsMotherEta(indexTrack);
593 if((fMotherEta != -1) && ((TMath::Abs(fPdg)) == 11)) fIsFrom = kElectronFromEta;
594
d2af20c5 595 fIndexTrack = indexTrack;
596
02524e30 597 }
598
d2af20c5 599 // MC PID for tagged
600 if(fgUseMCPID) {
601 if(TMath::Abs(fPdg) != 11) return;
602 }
603
02524e30 604 // Look at trackPart
d2af20c5 605 fIsPartner = kFALSE;
02524e30 606 Int_t pdgPart = GetPdg(indexTrackPart);
607 if(TMath::Abs(pdgPart) == 11) {
608 Int_t labMotherPart = GetLabMother(indexTrackPart);
d2af20c5 609 if((labMotherPart == fLabMother) && (indexTrack != indexTrackPart) && (TMath::Abs(fPdg) == 11) && (fPdg*pdgPart < 0) && (fLabMother >=0) && (fLabMother < (((AliStack *)fMCEvent->Stack())->GetNtrack()))) fIsPartner = kTRUE;
610 // special case of c and b
611 Int_t motherCPart = IsMotherC(indexTrackPart);
612 if((motherCPart != -1) && (fIsFrom == kElectronFromC) && (fPdg*pdgPart < 0)){
613 fIsPartner = kTRUE;
614 }
615 Int_t motherBPart = IsMotherB(indexTrackPart);
616 if((motherBPart != -1) && (fIsFrom == kElectronFromB) && (fPdg*pdgPart < 0)){
617 fIsPartner = kTRUE;
618 }
02524e30 619 }
d2af20c5 620
621 // Look at splitted tracks
622 fIsSplittedTrack = kFALSE;
623 if(indexTrackPart == fIndexTrack) fIsSplittedTrack = kTRUE;
624
02524e30 625 }
70da6c5a 626
02524e30 627 //////////////////////
628 // Sign
629 /////////////////////
630 Int_t sign = -1;
631 if((track->Charge() > 0.0) && (trackPart->Charge() > 0.0)) sign = kPp;
632 if((track->Charge() < 0.0) && (trackPart->Charge() < 0.0)) sign = kNn;
633 if(((track->Charge() > 0.0) && (trackPart->Charge() < 0.0)) || ((track->Charge() < 0.0) && (trackPart->Charge() > 0.0))) sign = kOs;
70da6c5a 634
635 /////////////////////////
636 // Cut effects
637 ////////////////////////
638 Double_t cuteffect[3];
639
640 if(fDebugLevel > 0) {
641 if(HasMCData()) {
642 cuteffect[0] = fPtESD;
643 cuteffect[1] = 0.0;
644 cuteffect[2] = fIsFrom;
69ac0e6f 645 if(!fCuts[0]){
646 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCe)))) fhtmp->Fill(cuteffect);
647 //if(fList->At(kMCe)) (dynamic_cast<THnSparseF *>(fList->At(kMCe)))->Fill(cuteffect);
70da6c5a 648 fCuts[0] = kTRUE;
649 }
650 }
651 }
652
653
654 ///////////////////////////////
655 // Cut effect: Partner track
656 ///////////////////////////////
657
658 if(fDebugLevel > 0) {
659 if(HasMCData() && fIsPartner) {
660 cuteffect[1] = 1.0;
661 if(!fCuts[1]) {
69ac0e6f 662 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCe)))) fhtmp->Fill(cuteffect);
663 //if(fList->At(kMCe)) (dynamic_cast<THnSparseF *>(fList->At(kMCe)))->Fill(cuteffect);
70da6c5a 664 fCuts[1] = kTRUE;
665 }
666 }
667 }
668
669 ///////////////////////////////
670 // Cut effect: Opposite sign
671 ///////////////////////////////
672
673 if(fDebugLevel > 0) {
674 if(HasMCData() && fIsPartner && (sign == kOs)) {
675 cuteffect[1] = 2.0;
676 if(!fCuts[2]) {
69ac0e6f 677 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCe)))) fhtmp->Fill(cuteffect);
678 //if(fList->At(kMCe)) (dynamic_cast<THnSparseF *>(fList->At(kMCe)))->Fill(cuteffect);
70da6c5a 679 fCuts[2] = kTRUE;
680 }
681 }
682 }
683
684 ////////////////////////
685 // Partner track cut
686 ////////////////////////
687 if(!SingleTrackCut(trackPart)) return;
688
689 if(fDebugLevel > 0) {
690 if(HasMCData() && fIsPartner && (sign==kOs)) {
691 cuteffect[1] = 3.0;
692 if(!fCuts[3]) {
69ac0e6f 693 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCe)))) fhtmp->Fill(cuteffect);
694 //if(fList->At(kMCe)) (dynamic_cast<THnSparseF *>(fList->At(kMCe)))->Fill(cuteffect);
70da6c5a 695 fCuts[3] = kTRUE;
696 }
697 }
698 }
699
700 /////////////////////////
701 // shared clusters cut
702 /////////////////////////
703 if(fSharedClusterCut && ShareCluster(track,trackPart)) return;
704
705 if(fDebugLevel > 0) {
706 if(HasMCData() && fIsPartner && (sign==kOs)) {
707 cuteffect[1] = 4.0;
708 if(!fCuts[4]) {
69ac0e6f 709 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCe)))) fhtmp->Fill(cuteffect);
710 //if(fList->At(kMCe)) (dynamic_cast<THnSparseF *>(fList->At(kMCe)))->Fill(cuteffect);
70da6c5a 711 fCuts[4] = kTRUE;
712 }
713 }
714 }
715
716 ////////////////////////
717 // PID Partner track
718 ////////////////////////
719 if(!PIDTrackCut(trackPart)) return;
720
721 if(fDebugLevel > 0) {
722 if(HasMCData() && fIsPartner && (sign==kOs)) {
723 cuteffect[1] = 5.0;
724 if(!fCuts[5]) {
69ac0e6f 725 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCe)))) fhtmp->Fill(cuteffect);
726 //if(fList->At(kMCe)) (dynamic_cast<THnSparseF *>(fList->At(kMCe)))->Fill(cuteffect);
70da6c5a 727 fCuts[5] = kTRUE;
728 }
729 }
730 }
731
732
02524e30 733 //////////////////////
734 // DCA
735 /////////////////////
736
737 Double_t xthis,xp;
738 Double_t dca = track->GetDCA(trackPart,fBz,xthis,xp);
69ac0e6f 739 if((fhtmpp = dynamic_cast<TH1F *>(fList->At(kDatadca)))) fhtmpp->Fill(dca);
70da6c5a 740 if(HasMCData()) {
741 //printf("has MC data for DCA\n");
742 //printf("IsPartner %d and isfrom %d\n",fIsPartner,fIsFrom);
69ac0e6f 743 if(fIsFrom==kElectronFromBackground) {
744 if((fhtmpf = dynamic_cast<TH2F *>(fList->At(kMCdca)))) fhtmpf->Fill(dca,fIsFrom);
745 }
70da6c5a 746 else {
69ac0e6f 747 if(fIsPartner){
748 if((fhtmpf = dynamic_cast<TH2F *>(fList->At(kMCdca)))) fhtmpf->Fill(dca,fIsFrom);
749 }
70da6c5a 750 }
751 }
69ac0e6f 752
02524e30 753 if(TMath::Abs(dca) > 3.0) return;
69ac0e6f 754
70da6c5a 755 if(fDebugLevel > 0) {
756 if(HasMCData() && fIsPartner && (sign==kOs)) {
757 cuteffect[1] = 6.0;
758 if(!fCuts[6]) {
69ac0e6f 759 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCe)))) fhtmp->Fill(cuteffect);
760 //if(fList->At(kMCe)) (dynamic_cast<THnSparseF *>(fList->At(kMCe)))->Fill(cuteffect);
70da6c5a 761 fCuts[6] = kTRUE;
762 }
763 }
764 }
765
02524e30 766 ///////////////////////////////////
767 // Calcul mother variables
768 ///////////////////////////////////
769 Double_t results[5];
770 Double_t resultsr[5];
70da6c5a 771
772
773 if(!fUseAliKFCode) {
774
775 /////////////////////////////
776 // Propagate
777 ////////////////////////////
778
779 Double_t norradius = TMath::Sqrt(fkVertex->GetX()*fkVertex->GetX()+fkVertex->GetY()*fkVertex->GetY());
780
bf892a6a 781 AliESDtrack trackCopy = AliESDtrack(*track);
782 AliESDtrack trackPartCopy = AliESDtrack(*trackPart);
70da6c5a 783 Bool_t propagateok = kTRUE;
bf892a6a 784 if((!(trackPartCopy.PropagateTo(norradius,fBz))) || (!(trackCopy.PropagateTo(norradius,fBz)))) propagateok = kFALSE;
70da6c5a 785 if(!propagateok) {
bf892a6a 786 //if(trackCopy) delete trackCopy;
787 //if(trackPartCopy) delete trackPartCopy;
70da6c5a 788 return;
789 }
790
bf892a6a 791 CalculateMotherVariable(&trackCopy,&trackPartCopy,&results[0]);
792 CalculateMotherVariableR(&trackCopy,&trackPartCopy,&resultsr[0]);
70da6c5a 793
bf892a6a 794 //if(trackCopy) delete trackCopy;
795 //if(trackPartCopy) delete trackPartCopy;
70da6c5a 796
797 }
798 else {
799
800 if(!CalculateMotherVariable(track,trackPart,&results[0])) return;
801 if(fDebugLevel > 0) {
802 if(HasMCData() && fIsPartner && (sign==kOs)) {
803 cuteffect[1] = 7.0;
804 if(!fCuts[7]) {
69ac0e6f 805 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCe)))) fhtmp->Fill(cuteffect);
806 //if(fList->At(kMCe)) (dynamic_cast<THnSparseF *>(fList->At(kMCe)))->Fill(cuteffect);
70da6c5a 807 fCuts[7] = kTRUE;
808 }
809 }
810 }
02524e30 811
70da6c5a 812 }
02524e30 813
814 /////////////////////////////////////
815 // Fill
816 /////////////////////////////////////
817
818 FillOutput(results, resultsr, sign);
70da6c5a 819
820 if(fDebugLevel > 0) {
821 if(HasMCData() && fIsPartner && (sign==kOs)) {
822
823 if(TMath::Abs(results[4]) < fOpeningAngleCut) {
824
825 cuteffect[1] = 8.0;
826 if(!fCuts[8]) {
69ac0e6f 827 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCe)))) fhtmp->Fill(cuteffect);
828 //if(fList->At(kMCe)) (dynamic_cast<THnSparseF *>(fList->At(kMCe)))->Fill(cuteffect);
70da6c5a 829 fCuts[8] = kTRUE;
830 }
831 if(TMath::Abs(results[1]) < fInvMassCut) {
832 cuteffect[1] = 9.0;
833 if(!fCuts[9]) {
69ac0e6f 834 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCe)))) fhtmp->Fill(cuteffect);
835 //if(fList->At(kMCe)) (dynamic_cast<THnSparseF *>(fList->At(kMCe)))->Fill(cuteffect);
70da6c5a 836 fCuts[9] = kTRUE;
837 }
838 }
839 }
840 }
841 }
02524e30 842
70da6c5a 843
02524e30 844}
845//_____________________________________________________________________________________
70da6c5a 846Bool_t AliHFEelecbackground::CalculateMotherVariable(AliESDtrack* const track, AliESDtrack* const trackpart, Double_t *results)
02524e30 847{
848 //
849 // variables mother and take the pt of the track
850 //
851 // results contain: ptmother, invmass, etamother, phimother, openingangle
852 //
70da6c5a 853 // with a chi2Ndf cut for AliKF code
854 //
02524e30 855
70da6c5a 856 if(!fUseAliKFCode) {
857
858 TVector3 v3Dtagged;
859 TVector3 v3Dpart;
860
67fe7bd0 861 Double_t pxyz[3];
862 track->PxPyPz(&pxyz[0]);
70da6c5a 863 v3Dtagged.SetXYZ(pxyz[0],pxyz[1],pxyz[2]);
864
67fe7bd0 865 Double_t pxyzpart[3];
866 trackpart->PxPyPz(&pxyzpart[0]);
70da6c5a 867 v3Dpart.SetXYZ(pxyzpart[0],pxyzpart[1],pxyzpart[2]);
868
869
870 TVector3 motherrec = v3Dtagged + v3Dpart;
871
872 Double_t etaESDmother = motherrec.Eta();
873 Double_t ptESDmother = motherrec.Pt();
874 Double_t phiESDmother = motherrec.Phi();
875 if(phiESDmother > TMath::Pi()) phiESDmother = phiESDmother - (2*TMath::Pi());
876
877
878 // openinganglepropagated
879 Double_t openingangle = v3Dtagged.Angle(v3Dpart);
880
881 // invmass
882 Double_t pESD = TMath::Sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]+pxyz[2]*pxyz[2]);
883 Double_t pESDpart = TMath::Sqrt(pxyzpart[0]*pxyzpart[0]+pxyzpart[1]*pxyzpart[1]+pxyzpart[2]*pxyzpart[2]);
884
885 // e propagate
886 Double_t eESD = TMath::Sqrt(pESD*pESD+fgkMe*fgkMe);
887 Double_t eESDpart = TMath::Sqrt(pESDpart*pESDpart+fgkMe*fgkMe);
888
889 Double_t invmass = TMath::Sqrt((eESD+eESDpart)*(eESD+eESDpart)-(motherrec.Px()*motherrec.Px()+motherrec.Py()*motherrec.Py()+motherrec.Pz()*motherrec.Pz()));
890
891 if(!results) return kFALSE;
02524e30 892
70da6c5a 893 results[0] = ptESDmother;
894 results[1] = etaESDmother;
895 results[2] = phiESDmother;
896 results[3] = invmass;
897 results[4] = openingangle;
898
899 return kTRUE;
02524e30 900
02524e30 901 }
70da6c5a 902 else {
903
904 AliKFParticle pair;
905 pair.Initialize();
906
907 // pid
908 Int_t pid1 = -11;
909 if(track->Charge() > 0.0) pid1 = 11;
910 Int_t pid2 = -11;
911 if(trackpart->Charge() > 0.0) pid2 = 11;
912
913
914 // daughters
915 AliKFParticle kf(*track,pid1);
916 AliKFParticle kfpart(*trackpart,pid2);
917
918 pair.AddDaughter(kf);
919 pair.AddDaughter(kfpart);
920
921 // variables
922 Double_t openingangle = kf.GetAngle(kfpart);
923 Double_t chi2ndf = pair.GetChi2()/pair.GetNDF();
924 //Double_t decayLength = pair.GetDecayLength();
925 Double_t radius = pair.GetR();
926 //Double_t masserror = pair.GetErrMass()>0?pair.GetErrMass()/pair.GetMass():1000000;
927 Double_t ptpair = pair.GetPt();
928 Double_t etapair = pair.GetEta();
929 Double_t phipair = pair.GetPhi();
930 Double_t masspair = pair.GetMass();
931
932 // Put them
933 if(!results) return kFALSE;
934
935 results[0] = ptpair;
936 results[1] = etapair;
937 results[2] = phipair;
938 results[3] = masspair;
939 results[4] = openingangle;
940
941 // chi2Ndf cut
69ac0e6f 942 if((fhtmpp = dynamic_cast<TH1F *>(fList->At(kDatachi2Ndf)))) fhtmpp->Fill(chi2ndf);
943 //if(fList->At(kDatachi2Ndf)) (dynamic_cast<TH1F *>(fList->At(kDatachi2Ndf)))->Fill(chi2ndf);
70da6c5a 944 if(HasMCData()){
69ac0e6f 945 if(fIsFrom==kElectronFromBackground) {
946 if((fhtmpf = dynamic_cast<TH2F *>(fList->At(kMCchi2Ndf)))) fhtmpf->Fill(chi2ndf,fIsFrom);
947 }
70da6c5a 948 else {
69ac0e6f 949 if(fIsPartner){
950 if((fhtmpf = dynamic_cast<TH2F *>(fList->At(kMCchi2Ndf)))) fhtmpf->Fill(chi2ndf,fIsFrom);
951 }
70da6c5a 952 }
953 }
954 if(chi2ndf > fChi2NdfCut) return kFALSE;
955 else {
69ac0e6f 956 if((fhtmpp = dynamic_cast<TH1F *>(fList->At(kDatar)))) fhtmpp->Fill(radius);
957 //if(fList->At(kDatar)) (dynamic_cast<TH1F *>(fList->At(kDatar)))->Fill(radius);
70da6c5a 958 if(HasMCData()) {
69ac0e6f 959 if(fIsFrom==kElectronFromBackground) {
960 if((fhtmpf = dynamic_cast<TH2F *>(fList->At(kMCr)))) fhtmpf->Fill(radius,fIsFrom);
961 //if(fList->At(kMCr))) (dynamic_cast<TH2F *>(fList->At(kMCr)))->Fill(radius,fIsFrom);
962 }
70da6c5a 963 else {
69ac0e6f 964 if(fIsPartner) {
965 if((fhtmpf = dynamic_cast<TH2F *>(fList->At(kMCr)))) fhtmpf->Fill(radius,fIsFrom);
966 }
70da6c5a 967 }
968 }
969 return kTRUE;
970 }
971
972 }
02524e30 973
974}
975//_____________________________________________________________________________________
976void AliHFEelecbackground::CalculateMotherVariableR(AliESDtrack* const track, AliESDtrack* const trackpart, Double_t *results)
977{
978 //
979 // variables mother
980 //
981 // results contain: ptmother, invmass, etamother, phimother, openingangle
70da6c5a 982 // Implemented only for no AliKF
02524e30 983 //
02524e30 984
70da6c5a 985 if(!fUseAliKFCode) {
986
987 TVector3 v3Dtagged;
988 TVector3 v3Dpart;
989
67fe7bd0 990 Double_t pxyz[3];
991 track->PxPyPz(&pxyz[0]);
70da6c5a 992 v3Dtagged.SetXYZ(pxyz[0],pxyz[1],pxyz[2]);
67fe7bd0 993 Double_t pxyzpart[3];
994 trackpart->PxPyPz(&pxyzpart[0]);
70da6c5a 995 v3Dpart.SetXYZ(pxyzpart[0],pxyzpart[1],pxyzpart[2]);
996
997 // rotate the partner
998 v3Dpart.RotateZ(TMath::Pi());
999 v3Dpart.GetXYZ(pxyzpart);
1000
1001
1002 TVector3 motherrec = v3Dtagged + v3Dpart;
1003
1004 Double_t etaESDmother = motherrec.Eta();
1005 Double_t ptESDmother = motherrec.Pt();
1006 Double_t phiESDmother = motherrec.Phi();
1007 if(phiESDmother > TMath::Pi()) phiESDmother = phiESDmother - (2*TMath::Pi());
1008
1009
1010 // openinganglepropagated
1011 Double_t openingangle = v3Dtagged.Angle(v3Dpart);
1012 //printf("Openingangle %f\n",openingangle);
1013
1014 // invmass
1015 Double_t pESD = TMath::Sqrt(pxyz[0]*pxyz[0]+pxyz[1]*pxyz[1]+pxyz[2]*pxyz[2]);
1016 Double_t pESDpart = TMath::Sqrt(pxyzpart[0]*pxyzpart[0]+pxyzpart[1]*pxyzpart[1]+pxyzpart[2]*pxyzpart[2]);
1017 // e propagate
1018 Double_t eESD = TMath::Sqrt(pESD*pESD+fgkMe*fgkMe);
1019 Double_t eESDpart = TMath::Sqrt(pESDpart*pESDpart+fgkMe*fgkMe);
1020
1021 Double_t invmass = TMath::Sqrt((eESD+eESDpart)*(eESD+eESDpart)-(motherrec.Px()*motherrec.Px()+motherrec.Py()*motherrec.Py()+motherrec.Pz()*motherrec.Pz()));
1022
1023 if(!results) return;
1024
1025 results[0] = ptESDmother;
1026 results[1] = etaESDmother;
1027 results[2] = phiESDmother;
1028 results[3] = invmass;
1029 results[4] = openingangle;
1030
02524e30 1031 }
02524e30 1032
1033}
1034//_________________________________________________________________________________
1035void AliHFEelecbackground::FillOutput(Double_t *results, Double_t *resultsr, Int_t sign)
1036{
1037 //
70da6c5a 1038 // Fill the Data and MC THnSparseF
02524e30 1039 //
70da6c5a 1040
1041 if((!results) || (!resultsr)) return;
02524e30 1042
70da6c5a 1043 Double_t co[6];
02524e30 1044 co[0] = fPtESD;
70da6c5a 1045 co[1] = results[0];
1046 co[2] = TMath::Abs(results[4]);
1047 co[3] = results[3];
1048 co[4] = sign;
1049 co[5] = 0.0;
02524e30 1050
69ac0e6f 1051 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kDatai)))) fhtmp->Fill(co);
1052 //if(fList->At(kDatai))(dynamic_cast<THnSparseF *>(fList->At(kDatai)))->Fill(co);
02524e30 1053
70da6c5a 1054 if((sign==kOs) && (!fUseAliKFCode)) {
1055
1056 co[1] = resultsr[0];
1057 co[2] = TMath::Abs(resultsr[4]);
1058 co[3] = resultsr[3];
1059 co[4] = kR;
1060 co[5] = 0.0;
02524e30 1061
69ac0e6f 1062 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kDatai)))) fhtmp->Fill(co);
1063 //if(fList->At(kDatai)) (dynamic_cast<THnSparseF *>(fList->At(kDatai)))->Fill(co);
70da6c5a 1064
1065 }
1066
1067 if(HasMCData()){
02524e30 1068
70da6c5a 1069 // Reset
1070 co[1] = results[0];
1071 co[2] = TMath::Abs(results[4]);
1072 co[3] = results[3];
02524e30 1073
70da6c5a 1074 // Origin
1075 co[4] = kElectronFromBackground;
1076 if((sign==kOs) && fIsPartner) co[4] = fIsFrom;
02524e30 1077
70da6c5a 1078 // Splitted tracks
1079 co[5] = kNotSplitted;
1080 if(fIsSplittedTrack) {
1081 if(sign==kOs){
1082 co[5] = kSplittedOs;
1083 }
1084 else {
1085 co[5] = kSplittedSs;
1086 }
1087 }
1088
69ac0e6f 1089 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCo)))) fhtmp->Fill(co);
1090 //if(fList->At(kMCo)) (dynamic_cast<THnSparseF *>(fList->At(kMCo)))->Fill(co);
70da6c5a 1091
02524e30 1092 }
70da6c5a 1093
1094}
02524e30 1095//_______________________________________________________________________________________________
70da6c5a 1096Bool_t AliHFEelecbackground::SingleTrackCut(AliESDtrack* const trackPart) const
02524e30 1097{
1098 //
1099 // Return minimum quality for the partner
1100 //
1101
70da6c5a 1102 if(trackPart->GetKinkIndex(0)>0) return kFALSE;
d2af20c5 1103
70da6c5a 1104
1105 UInt_t status = trackPart->GetStatus();
1106
1107 if(fRequireITSStandalone > 0) {
1108
1109 /////////////////////
1110 // ITS Standalone
1111 ////////////////////
1112
1113 if(fRequireITSStandalone==1) {
faee3b18 1114 if(((status & AliESDtrack::kITSin) == 0 || (trackPart->IsPureITSStandalone()) || ((status&AliESDtrack::kITSrefit)==0))) return kFALSE;
70da6c5a 1115 }
1116
1117 if(fRequireITSStandalone==2) {
1118 if(!trackPart->IsPureITSStandalone() || ((status&AliESDtrack::kITSrefit)==0)) return kFALSE;
1119 }
1120
1121 // Chi2
1122 Double_t chi2 = trackPart->GetITSchi2();
1123 if(chi2 > fMinITSChi2) return kFALSE;
1124
1125 // Min Nb of clusters
1126 Int_t nbcl = trackPart->GetITSclusters(0);
1127 if(nbcl < fMinNbCls) return kFALSE;
1128
1129 // Min Nb of points in SDD and SPD
1130 Int_t nbSDDSPD = 0;
1131 for(Int_t layer = 0; layer < 4; layer++){
1132 if(trackPart->HasPointOnITSLayer(layer)) nbSDDSPD++;
1133 }
1134 if(nbSDDSPD < fMinNbClsSDDSPD) return kFALSE;
1135
02524e30 1136
02524e30 1137 }
70da6c5a 1138 else {
1139
1140 /////////
1141 // TPC
1142 /////////
1143
1144 if((status&AliESDtrack::kTPCrefit)==0) return kFALSE;
1145
1146 // Min Nb of clusters
1147 Int_t nbcl = trackPart->GetTPCclusters(0);
1148 if(nbcl < fMinNbCls) return kFALSE;
1149
1150 }
1151
1152 return kTRUE;
1153
1154}
1155//_______________________________________________________________________________________________
1156Bool_t AliHFEelecbackground::PIDTrackCut(AliESDtrack* const trackPart)
1157{
1158 //
1159 // PID for the partner using TPC or ITS
1160 //
02524e30 1161
70da6c5a 1162 if(fRequireITSStandalone > 0) {
1163
1164 /////////////////////
1165 // ITS Standalone
1166 ////////////////////
1167
1168 // signal
1169 Double_t itsSignal = trackPart->GetITSsignal();
1170 Double_t p = trackPart->P();
1171
69ac0e6f 1172 if(fDebugLevel > 1) {
1173 if((fhtmpf = dynamic_cast<TH2F *>(fList->At(kMCcutPart0)))) fhtmpf->Fill(p,itsSignal);
1174 //if(fList->At(kMCcutPart0)) (dynamic_cast<TH2F *>(fList->At(kMCcutPart0)))->Fill(p,itsSignal);
70da6c5a 1175 }
1176
1177 ///////////
1178 // PID
1179 //////////
1180 if(fPIDPartner) {
1181
1182 // Take signal trackPart
1183 Double_t dEdxSamplesPart[4];
1184 trackPart->GetITSdEdxSamples(dEdxSamplesPart);
1185
1186 // Cut at 2 sigma
1187 if(!fPIDMethodPartnerITS) fPIDMethodPartnerITS = new AliESDpid;
1188 Float_t nsigma = fPIDMethodPartnerITS->NumberOfSigmasITS(trackPart, AliPID::kElectron);
1189 if(TMath::Abs(nsigma) > 2.0) return kFALSE;
1190
1191 // fill signal
1192 if(fDebugLevel > 1) {
1193
1194 Double_t entries[5];
1195 entries[0] = p;
1196 entries[1] = dEdxSamplesPart[0];
1197 entries[2] = dEdxSamplesPart[1];
1198 entries[3] = dEdxSamplesPart[2];
1199 entries[4] = dEdxSamplesPart[3];
1200
69ac0e6f 1201 //if(fList->At(kMCcutPart1)) (dynamic_cast<TH2F *>(fList->At(kMCcutPart1)))->Fill(p,itsSignal);
1202 if((fhtmpf = dynamic_cast<TH2F *>(fList->At(kMCcutPart1)))) fhtmpf->Fill(p,itsSignal);
1203 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCcutPart2)))) fhtmp->Fill(entries);
1204 //if(fList->At(kMCcutPart2)) (dynamic_cast<THnSparseF *>(fList->At(kMCcutPart2)))->Fill(entries);
70da6c5a 1205
1206 }
1207
1208 }
1209
1210 }
1211 else {
1212
1213 /////////
1214 // TPC
1215 /////////
1216
1217 Double_t tpcSignal = trackPart->GetTPCsignal();
1218 Double_t p = trackPart->GetInnerParam() ? trackPart->GetInnerParam()->P() : trackPart->P();
1219
1220 if(fDebugLevel > 1) {
1221 //printf("tpcSignal %f\n",tpcSignal);
69ac0e6f 1222 //if(fList->At(kMCcutPart0)) (dynamic_cast<TH2F *>(fList->At(kMCcutPart0)))->Fill(p,tpcSignal);
1223 if((fhtmpf = dynamic_cast<TH2F *>(fList->At(kMCcutPart0)))) fhtmpf->Fill(p,tpcSignal);
70da6c5a 1224 }
1225
1226 // PID
1227 if(fPIDPartner) {
1228 if(!fPIDMethodPartner) return kFALSE;
1229 AliHFEpidObject hfetrack;
3a72645a 1230 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
1231 hfetrack.SetRecTrack(trackPart);
70da6c5a 1232 //if(HasMCData()) hfetrack.fMCtrack = mctrack;
1233 if(!fPIDMethodPartner->IsSelected(&hfetrack)) return kFALSE;
1234
1235 if(fDebugLevel > 1) {
69ac0e6f 1236 if((fhtmpf = dynamic_cast<TH2F *>(fList->At(kMCcutPart1)))) fhtmpf->Fill(p,tpcSignal);
1237 //if(fList->At(kMCcutPart1)) (dynamic_cast<TH2F *>(fList->At(kMCcutPart1)))->Fill(p,tpcSignal);
70da6c5a 1238 }
1239 }
1240
1241 }
1242
1243 return kTRUE;
1244
1245}
1246//__________________________________________________________________________________________
69ac0e6f 1247Bool_t AliHFEelecbackground::ShareCluster(AliESDtrack * const track1,AliESDtrack * const track2)
70da6c5a 1248{
1249 //
1250 // Look if the two tracks shared clusters in the TPC or in the ITS depending on the method
1251 //
1252 // For TPC:
1253 // hsfval: number of shared clusters
1254 // hsmval: quality of the tracks
1255 //
1256 // For ITS:
1257 // compare the dEdx in the ITS
1258 //
1259
1260 if(!fRequireITSStandalone) {
1261
1262 //////////
1263 // TPC
1264 //////////
1265
1266
1267 Int_t nh = 0;
1268 Int_t an = 0;
1269 Int_t ns = 0;
1270 Float_t hsmval = 0.0;
1271 Float_t hsfval = 0.0;
1272
1273 for(unsigned int imap=0;imap<track1->GetTPCClusterMap().GetNbits(); imap++) {
1274 if(track1->GetTPCClusterMap().TestBitNumber(imap) &&
1275 track2->GetTPCClusterMap().TestBitNumber(imap)) {
1276 // Do they share it ?
1277 if (track1->GetTPCSharedMap().TestBitNumber(imap) &&
1278 track2->GetTPCSharedMap().TestBitNumber(imap))
1279 {
1280 an++;
1281 nh+=2;
1282 ns+=2;
1283 }
1284 else {
1285 an--;
1286 nh+=2;
1287 }
1288 }
1289 else if (track1->GetTPCClusterMap().TestBitNumber(imap) ||
1290 track2->GetTPCClusterMap().TestBitNumber(imap)) {
1291 an++;
1292 nh++;
1293 }
1294 }
1295
1296 if (nh >0) {
1297 hsmval = an*1.0/nh;
1298 hsfval = ns*1.0/nh;
1299 }
1300
1301
1302 if((hsfval > 0.15) || (hsmval > -0.15)) return kTRUE; //they share cluster
1303 else return kFALSE;
1304
02524e30 1305
1306 }
70da6c5a 1307 else {
1308
1309
1310 //////////
1311 // ITS
1312 /////////
1313
1314 // Take signals
1315 Double_t dEdxSamples1[4];
1316 track1->GetITSdEdxSamples(dEdxSamples1);
1317 Double_t dEdxSamples2[4];
1318 track2->GetITSdEdxSamples(dEdxSamples2);
1319
1320 // If there are matching
1321 Int_t nbClusters = 0;
1322 Bool_t match[4] = {kTRUE,kTRUE,kTRUE,kTRUE};
1323 Double_t limit[4] = {1.5,1.5,1.5,1.5};
1324 for(Int_t layer = 0; layer < 4; layer++) {
1325 if(track1->HasPointOnITSLayer(layer+2) && track2->HasPointOnITSLayer(layer+2)) {
1326 if(TMath::Abs(dEdxSamples1[layer]-dEdxSamples2[layer])>limit[layer]) match[layer] = kFALSE;
1327 nbClusters++;
1328 }
1329 }
1330 //printf("nbClusters %d\n",nbClusters);
1331
1332 // fill signal
1333 if(fDebugLevel > 1) {
1334 Double_t entriesSplit[5];
1335 entriesSplit[0] = 0.0;
1336 if(fIsSplittedTrack) entriesSplit[0] = 1.0;
1337
1338 for(Int_t layer = 0; layer < 4; layer++) {
1339 if(track1->HasPointOnITSLayer(layer+2) && track2->HasPointOnITSLayer(layer+2)) {
1340 entriesSplit[layer+1] = dEdxSamples1[layer]-dEdxSamples2[layer];
1341 }
1342 else entriesSplit[layer+1] = -100.0;
1343 }
69ac0e6f 1344 if((fhtmp = dynamic_cast<THnSparseF *>(fList->At(kMCcutPart3)))) fhtmp->Fill(entriesSplit);
1345 //if(fList->At(kMCcutPart3)) (dynamic_cast<THnSparseF *>(fList->At(kMCcutPart3)))->Fill(entriesSplit);
70da6c5a 1346 }
1347
1348 // Return
1349 Int_t nbClustersNotClose = 0;
1350 for(Int_t layer = 0; layer < 4; layer++) {
1351 if(!match[layer]) nbClustersNotClose++;
1352 }
1353 if((nbClusters > 1) && (nbClustersNotClose > 0.75*nbClusters)) return kFALSE;
1354 else return kTRUE;
1355
1356 }
02524e30 1357
70da6c5a 1358}
1359//____________________________________________________________________________________________________________
1360void AliHFEelecbackground::SetPIDPartner() {
1361
1362 //
1363 // Init the stuff for PID on the partner track
1364 //
1365
1366 fPIDPartner = kTRUE;
1367
1368 if(fRequireITSStandalone == 0) {
1369
1370 if(!fPIDMethodPartner) {
1371 fPIDMethodPartner = new AliHFEpid();
3a72645a 1372 fPIDMethodPartner->AddDetector("TPC", 0);
1373 fPIDMethodPartner->InitializePID(); // 3 sigma cut in TPC
70da6c5a 1374 }
1375
1376 }
1377 else {
1378
1379 if(!fPIDMethodPartnerITS) fPIDMethodPartnerITS = new AliESDpid;
1380
1381 }
02524e30 1382
1383}
1384//______________________________________________________________________________________________
1385void AliHFEelecbackground::SetEvent(AliESDEvent* const ESD)
1386{
1387 //
1388 // Set the AliESD Event, the magnetic field and the primary vertex
1389 //
1390
1391 fESD1=ESD;
1392 fBz=fESD1->GetMagneticField();
1393 fkVertex = fESD1->GetPrimaryVertex();
1394
1395}
1396//____________________________________________________________________________________________________________
1397Int_t AliHFEelecbackground::IsMotherGamma(Int_t tr) {
1398
1399 //
1400 // Return the lab of gamma mother or -1 if not gamma
1401 //
1402
1403 AliStack* stack = fMCEvent->Stack();
1404 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
1405
1406 // Take mother
1407 TParticle * particle = stack->Particle(tr);
1408 if(!particle) return -1;
1409 Int_t imother = particle->GetFirstMother();
1410 if((imother < 0) || (imother >= stack->GetNtrack())) return -1;
1411 TParticle * mother = stack->Particle(imother);
1412 if(!mother) return -1;
1413
1414 // Check gamma
1415 Int_t pdg = mother->GetPdgCode();
1416 if(TMath::Abs(pdg) == 22) return imother;
1417 if(TMath::Abs(pdg) == 11) {
1418 return IsMotherGamma(imother);
1419 }
1420 return -1;
1421
1422}
1423//
1424//____________________________________________________________________________________________________________
1425Int_t AliHFEelecbackground::IsMotherPi0(Int_t tr) {
1426
1427 //
1428 // Return the lab of pi0 mother or -1 if not pi0
1429 //
1430
1431 AliStack* stack = fMCEvent->Stack();
1432 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
1433
1434 // Take mother
1435 TParticle * particle = stack->Particle(tr);
1436 if(!particle) return -1;
1437 Int_t imother = particle->GetFirstMother();
1438 if((imother < 0) || (imother >= stack->GetNtrack())) return -1;
1439 TParticle * mother = stack->Particle(imother);
1440 if(!mother) return -1;
1441
1442 // Check gamma
1443 Int_t pdg = mother->GetPdgCode();
1444 if(TMath::Abs(pdg) == 111) return imother;
1445 if(TMath::Abs(pdg) == 11) {
1446 return IsMotherPi0(imother);
1447 }
1448 return -1;
1449
1450}
1451//____________________________________________________________________________________________________________
1452Int_t AliHFEelecbackground::IsMotherEta(Int_t tr) {
1453
1454 //
1455 // Return the lab of pi0 mother or -1 if not pi0
1456 //
1457
1458 AliStack* stack = fMCEvent->Stack();
1459 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
1460
1461 // Take mother
1462 TParticle * particle = stack->Particle(tr);
1463 if(!particle) return -1;
1464 Int_t imother = particle->GetFirstMother();
1465 if((imother < 0) || (imother >= stack->GetNtrack())) return -1;
1466 TParticle * mother = stack->Particle(imother);
1467 if(!mother) return -1;
1468
1469 // Check gamma
1470 Int_t pdg = mother->GetPdgCode();
1471 if(TMath::Abs(pdg) == 221) return imother;
1472 if(TMath::Abs(pdg) == 11) {
1473 return IsMotherEta(imother);
1474 }
1475 return -1;
1476
1477}
1478//____________________________________________________________________________________________________________
1479Int_t AliHFEelecbackground::IsMotherC(Int_t tr) {
1480
1481 //
1482 // Return the lab of signal mother or -1 if not signal
1483 //
1484
1485 AliStack* stack = fMCEvent->Stack();
1486 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
1487
1488 // Take mother
1489 TParticle * particle = stack->Particle(tr);
1490 if(!particle) return -1;
1491 Int_t imother = particle->GetFirstMother();
1492 if((imother < 0) || (imother >= stack->GetNtrack())) return -1;
1493 TParticle * mother = stack->Particle(imother);
1494 if(!mother) return -1;
1495
1496 // Check gamma
1497 Int_t pdg = mother->GetPdgCode();
1498 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;
1499 if(TMath::Abs(pdg) == 11) {
1500 return IsMotherC(imother);
1501 }
1502 return -1;
1503
1504}
1505//____________________________________________________________________________________________________________
1506Int_t AliHFEelecbackground::IsMotherB(Int_t tr) {
1507
1508 //
1509 // Return the lab of signal mother or -1 if not signal
1510 //
1511
1512 AliStack* stack = fMCEvent->Stack();
1513 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
1514
1515 // Take mother
1516 TParticle * particle = stack->Particle(tr);
1517 if(!particle) return -1;
1518 Int_t imother = particle->GetFirstMother();
1519 if((imother < 0) || (imother >= stack->GetNtrack())) return -1;
1520 TParticle * mother = stack->Particle(imother);
1521 if(!mother) return -1;
1522
1523 // Check gamma
1524 Int_t pdg = mother->GetPdgCode();
1525 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;
1526 if(TMath::Abs(pdg) == 11) {
1527 return IsMotherB(imother);
1528 }
1529 return -1;
1530
1531}
1532//____________________________________________________________________________________________________________
1533Int_t AliHFEelecbackground::GetPdg(Int_t tr) {
1534
1535 //
1536 // Simply pdg code
1537 //
1538
1539 AliStack* stack = fMCEvent->Stack();
1540 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
1541
1542 // MC Information
1543 TParticle * particle = stack->Particle(tr);
1544 if(!particle) return -1;
1545 Int_t pdg = particle->GetPdgCode();
1546
1547 return pdg;
1548
1549}
1550//____________________________________________________________________________________________________________
1551Int_t AliHFEelecbackground::GetLabMother(Int_t tr) {
1552
1553 //
1554 // Simply lab mother
1555 //
1556
1557 AliStack* stack = fMCEvent->Stack();
1558 if((tr < 0) || (tr >= stack->GetNtrack())) return -1;
1559
1560 // MC Information
1561 TParticle * particle = stack->Particle(tr);
1562 if(!particle) return -1;
1563 Int_t imother = particle->GetFirstMother();
1564
1565 return imother;
1566
1567}
1568//_______________________________________________________________________________________________
1569void AliHFEelecbackground::PostProcess()
1570{
1571 //
1572 // Post process the histos and extract the background pt spectra
1573 //
1574
1575 if(!fList) return;
1576
70da6c5a 1577 gStyle->SetPalette(1);
1578 gStyle->SetOptStat(1111);
1579 gStyle->SetPadBorderMode(0);
1580 gStyle->SetCanvasColor(10);
1581 gStyle->SetPadLeftMargin(0.13);
1582 gStyle->SetPadRightMargin(0.13);
1583
1584 /////////////////////////
1585 // Take the THnSparseF
1586 /////////////////////////
1587 THnSparseF *hsSparseData = dynamic_cast<THnSparseF *>(fList->FindObject("OpeningangleinvmassData"));
1588 THnSparseF *hsSparseMC = dynamic_cast<THnSparseF *>(fList->FindObject("OpeningangleinvmassMC"));
1589 THnSparseF *hsSparseCutPassedMC = dynamic_cast<THnSparseF *>(fList->FindObject("CutPassedMC"));
1590
1591 /////////////////////////////////
1592 // Cuts on the opening angle
1593 ////////////////////////////////
69ac0e6f 1594 if(!hsSparseData) return;
70da6c5a 1595 TAxis *axisOpeningAngleData = hsSparseData->GetAxis(2);
1596 Int_t binCutData = axisOpeningAngleData->FindBin(fOpeningAngleCut);
1597 hsSparseData->GetAxis(2)->SetRange(1,binCutData);
1598
1599 if(hsSparseMC) {
1600 TAxis *axisOpeningAngleMC = hsSparseMC->GetAxis(2);
1601 Int_t binCutMC = axisOpeningAngleMC->FindBin(fOpeningAngleCut);
1602 hsSparseMC->GetAxis(2)->SetRange(1,binCutMC);
1603 }
1604
1605 /////////////////////////
1606 // Prepare the histos
1607 ////////////////////////
1608
1609 TAxis *ptaxisinvmass = hsSparseData->GetAxis(3);
02524e30 1610 Int_t nbinsptinvmass = ptaxisinvmass->GetNbins();
1611
02524e30 1612 TH1D **invmassosptproj = new TH1D*[nbinsptinvmass];
1613 TH1D **invmassssptproj = new TH1D*[nbinsptinvmass];
1614 TH1D **invmassrptproj = new TH1D*[nbinsptinvmass];
1615 TH1D **invmassdiffptproj = new TH1D*[nbinsptinvmass];
1616 TH1D **invmassgammaptproj = new TH1D*[nbinsptinvmass];
1617 TH1D **invmasspi0ptproj = new TH1D*[nbinsptinvmass];
1618 TH1D **invmassetaptproj = new TH1D*[nbinsptinvmass];
1619 TH1D **invmassCptproj = new TH1D*[nbinsptinvmass];
1620 TH1D **invmassBptproj = new TH1D*[nbinsptinvmass];
1621
70da6c5a 1622 TH1D *yieldPtFound = (TH1D *) hsSparseData->Projection(0);
02524e30 1623 yieldPtFound->SetName("Found yield");
1624 yieldPtFound->Reset();
1625
1626 TH1D *yieldPtSourcesMC = 0x0;
70da6c5a 1627 TH1D *yieldPtSignalCutMC = 0x0;
1628 if(hsSparseMC) {
1629 yieldPtSourcesMC = (TH1D *) hsSparseMC->Projection(0);
02524e30 1630 yieldPtSourcesMC->SetName("Found yield");
1631 yieldPtSourcesMC->Reset();
70da6c5a 1632
1633 yieldPtSignalCutMC = (TH1D *) hsSparseMC->Projection(0);
02524e30 1634 yieldPtSignalCutMC->SetName("Found yield");
1635 yieldPtSignalCutMC->Reset();
1636 }
70da6c5a 1637
1638 ////////////
02524e30 1639 // canvas
70da6c5a 1640 ///////////
02524e30 1641 Int_t nbrow = (Int_t) (nbinsptinvmass/5);
70da6c5a 1642 TString namecanvas("InvMassSpectra");
02524e30 1643 TCanvas * canvas =new TCanvas(namecanvas,namecanvas,800,800);
1644 canvas->Divide(5,nbrow+1);
1645
70da6c5a 1646 /////////////////////////////
1647 // Loop on pt bins
1648 /////////////////////////////
02524e30 1649
02524e30 1650 for(Int_t k=1; k <= nbinsptinvmass; k++){
1651
1652 Double_t lowedge = ptaxisinvmass->GetBinLowEdge(k);
1653 Double_t upedge = ptaxisinvmass->GetBinUpEdge(k);
70da6c5a 1654
1655 // Pt bin
1656 hsSparseData->GetAxis(0)->SetRange(k,k);
1657 if(hsSparseMC) hsSparseMC->GetAxis(0)->SetRange(k,k);
02524e30 1658
70da6c5a 1659 //
1660 hsSparseData->GetAxis(4)->SetRange(kOs+1,kOs+1);
1661 invmassosptproj[k-1] = hsSparseData->Projection(3);
1662 hsSparseData->GetAxis(4)->SetRange(kPp+1,kNn+1);
1663 invmassssptproj[k-1] = hsSparseData->Projection(3);
1664 hsSparseData->GetAxis(4)->SetRange(kR+1,kR+1);
1665 invmassrptproj[k-1] = hsSparseData->Projection(3);
1666 hsSparseData->GetAxis(4)->SetRange(1,hsSparseData->GetAxis(4)->GetNbins());
02524e30 1667 invmassgammaptproj[k-1] = 0x0;
1668 invmasspi0ptproj[k-1] = 0x0;
1669 invmassetaptproj[k-1] = 0x0;
1670 invmassCptproj[k-1] = 0x0;
1671 invmassBptproj[k-1] = 0x0;
70da6c5a 1672 if(hsSparseMC) {
1673 hsSparseMC->GetAxis(4)->SetRange(kElectronFromGamma+1,kElectronFromGamma+1);
1674 invmassgammaptproj[k-1] = hsSparseMC->Projection(3);
1675 hsSparseMC->GetAxis(4)->SetRange(kElectronFromPi0+1,kElectronFromPi0+1);
1676 invmasspi0ptproj[k-1] = hsSparseMC->Projection(3);
1677 hsSparseMC->GetAxis(4)->SetRange(kElectronFromEta+1,kElectronFromEta+1);
1678 invmassetaptproj[k-1] = hsSparseMC->Projection(3);
1679 hsSparseMC->GetAxis(4)->SetRange(kElectronFromC+1,kElectronFromC+1);
1680 invmassCptproj[k-1] = hsSparseMC->Projection(3);
1681 hsSparseMC->GetAxis(4)->SetRange(kElectronFromB+1,kElectronFromB+1);
1682 invmassBptproj[k-1] = hsSparseMC->Projection(3);
1683 hsSparseMC->GetAxis(4)->SetRange(1,hsSparseMC->GetAxis(4)->GetNbins());
1684 }
1685
02524e30 1686 invmassdiffptproj[k-1] = (TH1D *) invmassosptproj[k-1]->Clone();
1687 TString name("Invmassdiffptbin");
1688 name += k;
1689 invmassdiffptproj[k-1]->SetName(name);
1690 invmassdiffptproj[k-1]->Add(invmassssptproj[k-1],-1.0);
1691
1692 TString namee("p_{T} tagged from ");
1693 namee += lowedge;
1694 namee += " GeV/c to ";
1695 namee += upedge;
1696 namee += " GeV/c";
1697
1698 invmassosptproj[k-1]->SetTitle((const char*)namee);
1699 invmassssptproj[k-1]->SetTitle((const char*)namee);
1700 invmassrptproj[k-1]->SetTitle((const char*)namee);
1701 invmassdiffptproj[k-1]->SetTitle((const char*)namee);
1702 if(invmassgammaptproj[k-1]) invmassgammaptproj[k-1]->SetTitle((const char*)namee);
1703 if(invmasspi0ptproj[k-1]) invmasspi0ptproj[k-1]->SetTitle((const char*)namee);
1704 if(invmassetaptproj[k-1]) invmassetaptproj[k-1]->SetTitle((const char*)namee);
1705 if(invmassCptproj[k-1]) invmassCptproj[k-1]->SetTitle((const char*)namee);
1706 if(invmassBptproj[k-1]) invmassBptproj[k-1]->SetTitle((const char*)namee);
70da6c5a 1707
02524e30 1708 invmassosptproj[k-1]->SetStats(0);
1709 invmassssptproj[k-1]->SetStats(0);
1710 invmassrptproj[k-1]->SetStats(0);
1711 invmassdiffptproj[k-1]->SetStats(0);
1712 if(invmassgammaptproj[k-1]) invmassgammaptproj[k-1]->SetStats(0);
1713 if(invmasspi0ptproj[k-1]) invmasspi0ptproj[k-1]->SetStats(0);
1714 if(invmassetaptproj[k-1]) invmassetaptproj[k-1]->SetStats(0);
1715 if(invmassCptproj[k-1]) invmassCptproj[k-1]->SetStats(0);
1716 if(invmassBptproj[k-1]) invmassBptproj[k-1]->SetStats(0);
1717
1718 Double_t yieldf = invmassdiffptproj[k-1]->Integral();
1719 if(invmassetaptproj[k-1] && invmasspi0ptproj[k-1] && invmassgammaptproj[k-1] && invmassCptproj[k-1] && invmassBptproj[k-1]) {
1720 Double_t yieldg = invmassetaptproj[k-1]->Integral() + invmasspi0ptproj[k-1]->Integral() + invmassgammaptproj[k-1]->Integral();
69ac0e6f 1721 if(yieldPtSourcesMC) yieldPtSourcesMC->SetBinContent(k,yieldg);
02524e30 1722
1723 Double_t yieldsignal = invmassCptproj[k-1]->Integral() + invmassBptproj[k-1]->Integral();
69ac0e6f 1724 if(yieldPtSignalCutMC) yieldPtSignalCutMC->SetBinContent(k,yieldsignal);
02524e30 1725 }
1726
1727 yieldPtFound->SetBinContent(k,yieldf);
1728
1729 canvas->cd(k);
1730 invmassosptproj[k-1]->Draw();
1731 invmassssptproj[k-1]->Draw("same");
1732 invmassdiffptproj[k-1]->Draw("same");
1733 invmassrptproj[k-1]->Draw("same");
1734 TLegend *legiv = new TLegend(0.4,0.6,0.89,0.89);
1735 legiv->AddEntry(invmassosptproj[k-1],"Opposite signs","p");
1736 legiv->AddEntry(invmassssptproj[k-1],"Same signs","p");
1737 legiv->AddEntry(invmassdiffptproj[k-1],"(Opposite - Same) signs","p");
1738 legiv->AddEntry(invmassrptproj[k-1],"rotated","p");
1739 if(invmassgammaptproj[k-1]) legiv->AddEntry(invmassgammaptproj[k-1],"e^{+}e^{-} from #gamma","p");
1740 if(invmasspi0ptproj[k-1]) legiv->AddEntry(invmasspi0ptproj[k-1],"e^{+}e^{-} from #pi^{0}","p");
1741 if(invmassetaptproj[k-1]) legiv->AddEntry(invmassetaptproj[k-1],"e^{+}e^{-} from #eta","p");
1742 legiv->Draw("same");
1743
70da6c5a 1744 hsSparseData->GetAxis(0)->SetRange(1,hsSparseData->GetAxis(0)->GetNbins());
1745 if(hsSparseMC) hsSparseMC->GetAxis(0)->SetRange(1,hsSparseMC->GetAxis(0)->GetNbins());
1746
02524e30 1747 }
1748
70da6c5a 1749 ////////////////////////////////////////////////////
1750 // End of plotting: do subtraction of background
1751 ///////////////////////////////////////////////////
1752
02524e30 1753 yieldPtFound->SetStats(0);
1754 if(yieldPtSourcesMC) yieldPtSourcesMC->SetStats(0);
1755 if(yieldPtSignalCutMC) yieldPtSignalCutMC->SetStats(0);
1756
70da6c5a 1757 TCanvas * canvasfin =new TCanvas("ResultsElecBackGround","ResultsElecBackGround",800,800);
02524e30 1758 canvasfin->cd(1);
1759 yieldPtFound->Draw();
1760 if(yieldPtSourcesMC && yieldPtSignalCutMC) {
1761 yieldPtSourcesMC->Draw("same");
1762 yieldPtSignalCutMC->Draw("same");
1763 TLegend *lega = new TLegend(0.4,0.6,0.89,0.89);
70da6c5a 1764 lega->AddEntry(yieldPtFound,"Contributions found","l");
1765 lega->AddEntry(yieldPtSourcesMC,"Contributions of e^{+}e^{-} from #gamma, #pi^{0} and #eta","l");
1766 lega->AddEntry(yieldPtSignalCutMC,"Contributions of e^{+}e^{-} from C and B","l");
02524e30 1767 lega->Draw("same");
1768 }
1769
70da6c5a 1770 if(hsSparseCutPassedMC){
1771 hsSparseCutPassedMC->GetAxis(1)->SetRange(1,1);
1772 hsSparseCutPassedMC->GetAxis(2)->SetRange(1,4);
1773 TH1D *hsSparseCutPassedMCproj = hsSparseCutPassedMC->Projection(0);
02524e30 1774
70da6c5a 1775 TH1D *cYieldPtFound = (TH1D*)yieldPtFound->Clone("RatioEfficiency");
1776 if(hsSparseCutPassedMCproj->Integral() > 0.0) cYieldPtFound->Divide(hsSparseCutPassedMCproj);
1777
1778 TCanvas * canvasfratio =new TCanvas("RatioEfficiency","RatioEfficiency",800,800);
1779 canvasfratio->cd(1);
1780 cYieldPtFound->Draw();
1781 }
1782
1783 //////////////////////////
1784 // fListPostProcess
1785 /////////////////////////
02524e30 1786
1787 if(!fListPostProcess) fListPostProcess = new TList();
1788 fListPostProcess->SetName("ListPostProcess");
1789
1790 for(Int_t k=0; k < nbinsptinvmass; k++){
1791 fListPostProcess->AddAt(invmassosptproj[k],kOos+kNOutput*k);
1792 fListPostProcess->AddAt(invmassssptproj[k],kOss+kNOutput*k);
1793 fListPostProcess->AddAt(invmassrptproj[k],kOr+kNOutput*k);
1794 fListPostProcess->AddAt(invmassdiffptproj[k],kOdiff+kNOutput*k);
1795 if(invmassgammaptproj[k]) fListPostProcess->AddAt(invmassgammaptproj[k],kNOutput*nbinsptinvmass+kNMCInfo*k+kElectronFromGamma);
1796 if(invmasspi0ptproj[k]) fListPostProcess->AddAt(invmasspi0ptproj[k],kNOutput*nbinsptinvmass+kNMCInfo*k+kElectronFromPi0);
1797 if(invmassetaptproj[k]) fListPostProcess->AddAt(invmassetaptproj[k],kNOutput*nbinsptinvmass+kNMCInfo*k+kElectronFromEta);
1798 if(invmassCptproj[k]) fListPostProcess->AddAt(invmassCptproj[k],kNOutput*nbinsptinvmass+kNMCInfo*k+kElectronFromC);
1799 if(invmassBptproj[k]) fListPostProcess->AddAt(invmassBptproj[k],kNOutput*nbinsptinvmass+kNMCInfo*k+kElectronFromB);
1800 }
1801
1802 fListPostProcess->AddAt(yieldPtFound,kNOutput*nbinsptinvmass+kNMCInfo*nbinsptinvmass);
1803 if(yieldPtSourcesMC) fListPostProcess->AddAt(yieldPtSourcesMC,kNOutput*nbinsptinvmass+kNMCInfo*nbinsptinvmass+1);
1804 if(yieldPtSignalCutMC) fListPostProcess->AddAt(yieldPtSignalCutMC,kNOutput*nbinsptinvmass+kNMCInfo*nbinsptinvmass+2);
1805
1806 // delete dynamic array
1807 delete[] invmassosptproj;
1808 delete[] invmassssptproj;
1809 delete[] invmassrptproj;
1810 delete[] invmassdiffptproj;
1811 delete[] invmassgammaptproj;
1812 delete[] invmasspi0ptproj;
1813 delete[] invmassetaptproj;
1814 delete[] invmassCptproj;
1815 delete[] invmassBptproj;
1816
1817}
1818//_______________________________________________________________________________________________
1819void AliHFEelecbackground::Plot() const
1820{
1821 //
1822 // Plot the output
1823 //
1824
1825 if(!fList) return;
70da6c5a 1826
1827 gStyle->SetPalette(1);
1828 gStyle->SetOptStat(1111);
1829 gStyle->SetPadBorderMode(0);
1830 gStyle->SetCanvasColor(10);
1831 gStyle->SetPadLeftMargin(0.13);
1832 gStyle->SetPadRightMargin(0.13);
1833
1834
1835 /////////////////////////
1836 // Take the THnSparseF
1837 /////////////////////////
1838 THnSparseF *hsSparseData = dynamic_cast<THnSparseF *>(fList->FindObject("OpeningangleinvmassData"));
1839 THnSparseF *hsSparseMC = dynamic_cast<THnSparseF *>(fList->FindObject("OpeningangleinvmassMC"));
69ac0e6f 1840 if(!hsSparseData) return;
70da6c5a 1841
1842 ////////////////////
1843 // Opening angle
1844 ////////////////////
1845
1846 // Opening angle one direction
1847 hsSparseData->GetAxis(4)->SetRange(kPp+1,kPp+1);
1848 TH1D *openingangleppproj = hsSparseData->Projection(2);
1849 hsSparseData->GetAxis(4)->SetRange(kNn+1,kNn+1);
1850 TH1D *openinganglennproj = hsSparseData->Projection(2);
1851 hsSparseData->GetAxis(4)->SetRange(kPp+1,kNn+1);
1852 TH1D *openinganglessproj = hsSparseData->Projection(2);
1853 hsSparseData->GetAxis(4)->SetRange(kR+1,kR+1);
1854 TH1D *openinganglerproj = hsSparseData->Projection(2);
1855 hsSparseData->GetAxis(4)->SetRange(kOs+1,kOs+1);
1856 TH1D *openingangleosproj = hsSparseData->Projection(2);
1857 hsSparseData->GetAxis(4)->SetRange(1,hsSparseData->GetAxis(4)->GetNbins());
02524e30 1858
1859 TH1D *openinganglegammaproj = 0x0;
1860 TH1D *openinganglepi0proj = 0x0;
1861 TH1D *openingangleCproj = 0x0;
1862 TH1D *openingangleBproj = 0x0;
1863 TH1D *openingangleetaproj = 0x0;
d2af20c5 1864 TH1D *openingangleSplittedTrackssproj = 0x0;
1865 TH1D *openingangleSplittedTrackosproj = 0x0;
70da6c5a 1866 if(hsSparseMC) {
1867 hsSparseMC->GetAxis(4)->SetRange(kElectronFromGamma+1,kElectronFromGamma+1);
1868 openinganglegammaproj = hsSparseMC->Projection(2);
1869 hsSparseMC->GetAxis(4)->SetRange(kElectronFromPi0+1,kElectronFromPi0+1);
1870 openinganglepi0proj = hsSparseMC->Projection(2);
1871 hsSparseMC->GetAxis(4)->SetRange(kElectronFromEta+1,kElectronFromEta+1);
1872 openingangleetaproj = hsSparseMC->Projection(2);
1873 hsSparseMC->GetAxis(4)->SetRange(kElectronFromC+1,kElectronFromC+1);
1874 openingangleCproj = hsSparseMC->Projection(2);
1875 hsSparseMC->GetAxis(4)->SetRange(kElectronFromB+1,kElectronFromB+1);
1876 openingangleBproj = hsSparseMC->Projection(2);
1877 hsSparseMC->GetAxis(4)->SetRange(1,hsSparseMC->GetAxis(4)->GetNbins());
1878 hsSparseMC->GetAxis(5)->SetRange(kSplittedSs+1,kSplittedSs+1);
1879 openingangleSplittedTrackssproj = hsSparseMC->Projection(2);
1880 hsSparseMC->GetAxis(5)->SetRange(kSplittedOs+1,kSplittedOs+1);
1881 openingangleSplittedTrackosproj = hsSparseMC->Projection(2);
1882 hsSparseMC->GetAxis(5)->SetRange(1,hsSparseMC->GetAxis(5)->GetNbins());
1883 }
1884
1885 // Projection pt-opening angle
1886 hsSparseData->GetAxis(4)->SetRange(kPp+1,kPp+1);
1887 TH2D *openingangleppproj2D = hsSparseData->Projection(0,2);
1888 hsSparseData->GetAxis(4)->SetRange(kNn+1,kNn+1);
1889 TH2D *openinganglennproj2D = hsSparseData->Projection(0,2);
1890 hsSparseData->GetAxis(4)->SetRange(kPp+1,kNn+1);
1891 TH2D *openinganglessproj2D = hsSparseData->Projection(0,2);
1892 hsSparseData->GetAxis(4)->SetRange(kR+1,kR+1);
1893 TH2D *openinganglerproj2D = hsSparseData->Projection(0,2);
1894 hsSparseData->GetAxis(4)->SetRange(kOs+1,kOs+1);
1895 TH2D *openingangleosproj2D = hsSparseData->Projection(0,2);
1896 hsSparseData->GetAxis(4)->SetRange(1,hsSparseData->GetAxis(4)->GetNbins());
02524e30 1897
70da6c5a 1898 TH2D *openinganglegammaproj2D = 0x0;
1899 TH2D *openinganglepi0proj2D = 0x0;
1900 TH2D *openingangleCproj2D = 0x0;
1901 TH2D *openingangleBproj2D = 0x0;
1902 TH2D *openingangleetaproj2D = 0x0;
1903 TH2D *openingangleSplittedTrackssproj2D = 0x0;
1904 TH2D *openingangleSplittedTrackosproj2D = 0x0;
1905 if(hsSparseMC) {
1906 hsSparseMC->GetAxis(4)->SetRange(kElectronFromGamma+1,kElectronFromGamma+1);
1907 openinganglegammaproj2D = hsSparseMC->Projection(0,2);
1908 hsSparseMC->GetAxis(4)->SetRange(kElectronFromPi0+1,kElectronFromPi0+1);
1909 openinganglepi0proj2D = hsSparseMC->Projection(0,2);
1910 hsSparseMC->GetAxis(4)->SetRange(kElectronFromEta+1,kElectronFromEta+1);
1911 openingangleetaproj2D = hsSparseMC->Projection(0,2);
1912 hsSparseMC->GetAxis(4)->SetRange(kElectronFromC+1,kElectronFromC+1);
1913 openingangleCproj2D = hsSparseMC->Projection(0,2);
1914 hsSparseMC->GetAxis(4)->SetRange(kElectronFromB+1,kElectronFromB+1);
1915 openingangleBproj2D = hsSparseMC->Projection(0,2);
1916 hsSparseMC->GetAxis(4)->SetRange(1, hsSparseMC->GetAxis(4)->GetNbins());
1917 hsSparseMC->GetAxis(5)->SetRange(kSplittedSs+1,kSplittedSs+1);
1918 openingangleSplittedTrackssproj2D = hsSparseMC->Projection(0,2);
1919 hsSparseMC->GetAxis(5)->SetRange(kSplittedOs+1,kSplittedOs+1);
1920 openingangleSplittedTrackosproj2D = hsSparseMC->Projection(0,2);
1921 hsSparseMC->GetAxis(5)->SetRange(1,hsSparseMC->GetAxis(5)->GetNbins());
1922 }
02524e30 1923
70da6c5a 1924 openingangleppproj2D->SetStats(0);
1925 openinganglennproj2D->SetStats(0);
1926 openinganglessproj2D->SetStats(0);
1927 openinganglerproj2D->SetStats(0);
1928 openingangleosproj2D->SetStats(0);
1929 if(openinganglegammaproj2D) openinganglegammaproj2D->SetStats(0);
1930 if(openinganglepi0proj2D) openinganglepi0proj2D->SetStats(0);
1931 if(openingangleCproj2D) openingangleCproj2D->SetStats(0);
1932 if(openingangleBproj2D) openingangleBproj2D->SetStats(0);
1933 if(openingangleetaproj2D) openingangleetaproj2D->SetStats(0);
1934 if(openingangleSplittedTrackssproj2D) openingangleSplittedTrackssproj2D->SetStats(0);
1935 if(openingangleSplittedTrackosproj2D) openingangleSplittedTrackosproj2D->SetStats(0);
1936
1937 openingangleppproj2D->SetTitle("openingangleppproj2D");
1938 openinganglennproj2D->SetTitle("openinganglennproj2D");
1939 openinganglessproj2D->SetTitle("openinganglessproj2D");
1940 openinganglerproj2D->SetTitle("openinganglerproj2D");
1941 openingangleosproj2D->SetTitle("openingangleosproj2D");
1942 if(openinganglegammaproj2D) openinganglegammaproj2D->SetTitle("openinganglegammaproj2D");
1943 if(openinganglepi0proj2D) openinganglepi0proj2D->SetTitle("openinganglepi0proj2D");
1944 if(openingangleCproj2D) openingangleCproj2D->SetTitle("openingangleCproj2D");
1945 if(openingangleBproj2D) openingangleBproj2D->SetTitle("openingangleBproj2D");
1946 if(openingangleetaproj2D) openingangleetaproj2D->SetTitle("openingangleetaproj2D");
1947 if(openingangleSplittedTrackssproj2D) openingangleSplittedTrackssproj2D->SetTitle("openingangleSplittedTrackssproj2D");
1948 if(openingangleSplittedTrackosproj2D) openingangleSplittedTrackosproj2D->SetTitle("openingangleSplittedTrackosproj2D");
02524e30 1949
02524e30 1950 openingangleppproj->SetStats(0);
1951 openinganglennproj->SetStats(0);
1952 openinganglessproj->SetStats(0);
1953 openinganglerproj->SetStats(0);
1954 openingangleosproj->SetStats(0);
1955 if(openinganglegammaproj) openinganglegammaproj->SetStats(0);
1956 if(openinganglepi0proj) openinganglepi0proj->SetStats(0);
1957 if(openingangleCproj) openingangleCproj->SetStats(0);
1958 if(openingangleBproj) openingangleBproj->SetStats(0);
1959 if(openingangleetaproj) openingangleetaproj->SetStats(0);
d2af20c5 1960 if(openingangleSplittedTrackssproj) openingangleSplittedTrackssproj->SetStats(0);
1961 if(openingangleSplittedTrackosproj) openingangleSplittedTrackosproj->SetStats(0);
02524e30 1962
1963 openingangleppproj->SetTitle("");
1964 openinganglennproj->SetTitle("");
1965 openinganglessproj->SetTitle("");
1966 openinganglerproj->SetTitle("");
1967 openingangleosproj->SetTitle("");
1968 if(openinganglegammaproj) openinganglegammaproj->SetTitle("");
1969 if(openinganglepi0proj) openinganglepi0proj->SetTitle("");
1970 if(openingangleCproj) openingangleCproj->SetTitle("");
1971 if(openingangleBproj) openingangleBproj->SetTitle("");
1972 if(openingangleetaproj) openingangleetaproj->SetTitle("");
d2af20c5 1973 if(openingangleSplittedTrackssproj) openingangleSplittedTrackssproj->SetTitle("");
1974 if(openingangleSplittedTrackosproj) openingangleSplittedTrackosproj->SetTitle("");
02524e30 1975
70da6c5a 1976 ////////////////////////////
1977 // Invariant mass
1978 ///////////////////////////
1979
1980 // Cuts on the opening angle
1981 TAxis *axisOpeningAngleData = hsSparseData->GetAxis(2);
1982 Int_t binCutData = axisOpeningAngleData->FindBin(fOpeningAngleCut);
1983 hsSparseData->GetAxis(2)->SetRange(1,binCutData);
1984
1985 // Debug
1986 //printf("Get Bin low edge %f, Get Bin Up edge %f for hsSparseData\n",axisOpeningAngleData->GetBinLowEdge(binCutData),axisOpeningAngleData->GetBinUpEdge(binCutData));
1987
1988 // Invariant mass
1989 hsSparseData->GetAxis(4)->SetRange(kPp+1,kPp+1);
1990 TH1D *invmassppproj = hsSparseData->Projection(3);
1991 hsSparseData->GetAxis(4)->SetRange(kNn+1,kNn+1);
1992 TH1D *invmassnnproj = hsSparseData->Projection(3);
1993 hsSparseData->GetAxis(4)->SetRange(kPp+1,kNn+1);
1994 TH1D *invmassssproj = hsSparseData->Projection(3);
1995 hsSparseData->GetAxis(4)->SetRange(kR+1,kR+1);
1996 TH1D *invmassrproj = hsSparseData->Projection(3);
1997 hsSparseData->GetAxis(4)->SetRange(kOs+1,kOs+1);
1998 TH1D *invmassosproj = hsSparseData->Projection(3);
1999 hsSparseData->GetAxis(4)->SetRange(1,hsSparseData->GetAxis(4)->GetNbins());
2000
2001 TH1D *invmassgammaproj = 0x0;
2002 TH1D *invmasspi0proj = 0x0;
2003 TH1D *invmassCproj = 0x0;
2004 TH1D *invmassBproj = 0x0;
2005 TH1D *invmassetaproj = 0x0;
2006 TH1D *invmassSplittedTrackssproj = 0x0;
2007 TH1D *invmassSplittedTrackosproj = 0x0;
2008 if(hsSparseMC) {
2009 TAxis *axisOpeningAngleMC = hsSparseMC->GetAxis(2);
2010 Int_t binCutMC = axisOpeningAngleMC->FindBin(fOpeningAngleCut);
2011 hsSparseMC->GetAxis(2)->SetRange(1,binCutMC);
2012
2013 // Debug
2014 //printf("Get Bin low edge %f, Get Bin Up edge %f for hsSparseMC\n",axisOpeningAngleMC->GetBinLowEdge(binCutMC),axisOpeningAngleMC->GetBinUpEdge(binCutMC));
2015
2016 hsSparseMC->GetAxis(4)->SetRange(kElectronFromGamma+1,kElectronFromGamma+1);
2017 invmassgammaproj = hsSparseMC->Projection(3);
2018 hsSparseMC->GetAxis(4)->SetRange(kElectronFromPi0+1,kElectronFromPi0+1);
2019 invmasspi0proj = hsSparseMC->Projection(3);
2020 hsSparseMC->GetAxis(4)->SetRange(kElectronFromEta+1,kElectronFromEta+1);
2021 invmassetaproj = hsSparseMC->Projection(3);
2022 hsSparseMC->GetAxis(4)->SetRange(kElectronFromC+1,kElectronFromC+1);
2023 invmassCproj = hsSparseMC->Projection(3);
2024 hsSparseMC->GetAxis(4)->SetRange(kElectronFromB+1,kElectronFromB+1);
2025 invmassBproj = hsSparseMC->Projection(3);
2026 hsSparseMC->GetAxis(4)->SetRange(1,hsSparseMC->GetAxis(4)->GetNbins());
2027 hsSparseMC->GetAxis(5)->SetRange(kSplittedSs+1,kSplittedSs+1);
2028 invmassSplittedTrackssproj = hsSparseMC->Projection(3);
2029 hsSparseMC->GetAxis(5)->SetRange(kSplittedOs+1,kSplittedOs+1);
2030 invmassSplittedTrackosproj = hsSparseMC->Projection(3);
2031 hsSparseMC->GetAxis(5)->SetRange(1,hsSparseMC->GetAxis(5)->GetNbins());
2032 }
2033
2034 invmassppproj->SetStats(0);
2035 invmassnnproj->SetStats(0);
2036 invmassssproj->SetStats(0);
2037 invmassrproj->SetStats(0);
2038 invmassosproj->SetStats(0);
2039 if(invmassgammaproj) invmassgammaproj->SetStats(0);
2040 if(invmasspi0proj) invmasspi0proj->SetStats(0);
2041 if(invmassCproj) invmassCproj->SetStats(0);
2042 if(invmassBproj) invmassBproj->SetStats(0);
2043 if(invmassetaproj) invmassetaproj->SetStats(0);
2044 if(invmassSplittedTrackssproj) invmassSplittedTrackssproj->SetStats(0);
2045 if(invmassSplittedTrackosproj) invmassSplittedTrackosproj->SetStats(0);
2046
02524e30 2047 invmassppproj->SetTitle("");
2048 invmassnnproj->SetTitle("");
2049 invmassssproj->SetTitle("");
2050 invmassrproj->SetTitle("");
2051 invmassosproj->SetTitle("");
2052 if(invmassgammaproj) invmassgammaproj->SetTitle("");
2053 if(invmasspi0proj) invmasspi0proj->SetTitle("");
2054 if(invmassCproj) invmassCproj->SetTitle("");
2055 if(invmassBproj) invmassBproj->SetTitle("");
2056 if(invmassetaproj) invmassetaproj->SetTitle("");
d2af20c5 2057 if(invmassSplittedTrackssproj) invmassSplittedTrackssproj->SetTitle("");
2058 if(invmassSplittedTrackosproj) invmassSplittedTrackosproj->SetTitle("");
2059
70da6c5a 2060 // Projection pt-invariant mass angle
2061 hsSparseData->GetAxis(4)->SetRange(kPp+1,kPp+1);
2062 TH2D *invmassppproj2D = hsSparseData->Projection(0,3);
2063 hsSparseData->GetAxis(4)->SetRange(kNn+1,kNn+1);
2064 TH2D *invmassnnproj2D = hsSparseData->Projection(0,3);
2065 hsSparseData->GetAxis(4)->SetRange(kPp+1,kNn+1);
2066 TH2D *invmassssproj2D = hsSparseData->Projection(0,3);
2067 hsSparseData->GetAxis(4)->SetRange(kR+1,kR+1);
2068 TH2D *invmassrproj2D = hsSparseData->Projection(0,3);
2069 hsSparseData->GetAxis(4)->SetRange(kOs+1,kOs+1);
2070 TH2D *invmassosproj2D = hsSparseData->Projection(0,3);
2071 hsSparseData->GetAxis(4)->SetRange(1,hsSparseData->GetAxis(4)->GetNbins());
d2af20c5 2072
2073 TH2D *invmassgammaproj2D = 0x0;
2074 TH2D *invmasspi0proj2D = 0x0;
2075 TH2D *invmassCproj2D = 0x0;
2076 TH2D *invmassBproj2D = 0x0;
2077 TH2D *invmassetaproj2D = 0x0;
2078 TH2D *invmassSplittedTrackssproj2D = 0x0;
2079 TH2D *invmassSplittedTrackosproj2D = 0x0;
70da6c5a 2080 if(hsSparseMC) {
2081 hsSparseMC->GetAxis(4)->SetRange(kElectronFromGamma+1,kElectronFromGamma+1);
2082 invmassgammaproj2D = hsSparseMC->Projection(0,3);
2083 hsSparseMC->GetAxis(4)->SetRange(kElectronFromPi0+1,kElectronFromPi0+1);
2084 invmasspi0proj2D = hsSparseMC->Projection(0,3);
2085 hsSparseMC->GetAxis(4)->SetRange(kElectronFromEta+1,kElectronFromEta+1);
2086 invmassetaproj2D = hsSparseMC->Projection(0,3);
2087 hsSparseMC->GetAxis(4)->SetRange(kElectronFromC+1,kElectronFromC+1);
2088 invmassCproj2D = hsSparseMC->Projection(0,3);
2089 hsSparseMC->GetAxis(4)->SetRange(kElectronFromB+1,kElectronFromB+1);
2090 invmassBproj2D = hsSparseMC->Projection(0,3);
2091 hsSparseMC->GetAxis(4)->SetRange(1,hsSparseMC->GetAxis(4)->GetNbins());
2092 hsSparseMC->GetAxis(5)->SetRange(kSplittedSs+1,kSplittedSs+1);
2093 invmassSplittedTrackssproj2D = hsSparseMC->Projection(0,3);
2094 hsSparseMC->GetAxis(5)->SetRange(kSplittedOs+1,kSplittedOs+1);
2095 invmassSplittedTrackosproj2D = hsSparseMC->Projection(0,3);
2096 hsSparseMC->GetAxis(5)->SetRange(1,hsSparseMC->GetAxis(5)->GetNbins());
2097 }
2098
d2af20c5 2099
2100 invmassppproj2D->SetStats(0);
2101 invmassnnproj2D->SetStats(0);
2102 invmassssproj2D->SetStats(0);
2103 invmassrproj2D->SetStats(0);
2104 invmassosproj2D->SetStats(0);
2105 if(invmassgammaproj2D) invmassgammaproj2D->SetStats(0);
2106 if(invmasspi0proj2D) invmasspi0proj2D->SetStats(0);
2107 if(invmassCproj2D) invmassCproj2D->SetStats(0);
2108 if(invmassBproj2D) invmassBproj2D->SetStats(0);
2109 if(invmassetaproj2D) invmassetaproj2D->SetStats(0);
2110 if(invmassSplittedTrackssproj2D) invmassSplittedTrackssproj2D->SetStats(0);
2111 if(invmassSplittedTrackosproj2D) invmassSplittedTrackosproj2D->SetStats(0);
70da6c5a 2112
d2af20c5 2113 invmassppproj2D->SetTitle("invmassppproj2D");
2114 invmassnnproj2D->SetTitle("invmassnnproj2D");
2115 invmassssproj2D->SetTitle("invmassssproj2D");
2116 invmassrproj2D->SetTitle("invmassrproj2D");
2117 invmassosproj2D->SetTitle("invmassosproj2D");
2118 if(invmassgammaproj2D) invmassgammaproj2D->SetTitle("invmassgammaproj2D");
2119 if(invmasspi0proj2D) invmasspi0proj2D->SetTitle("invmasspi0proj2D");
2120 if(invmassCproj2D) invmassCproj2D->SetTitle("invmassCproj2D");
2121 if(invmassBproj2D) invmassBproj2D->SetTitle("invmassBproj2D");
2122 if(invmassetaproj2D) invmassetaproj2D->SetTitle("invmassetaproj2D");
2123 if(invmassSplittedTrackssproj2D) invmassSplittedTrackssproj2D->SetTitle("invmassSplittedTrackssproj2D");
2124 if(invmassSplittedTrackosproj2D) invmassSplittedTrackosproj2D->SetTitle("invmassSplittedTrackosproj2D");
2125
02524e30 2126
70da6c5a 2127 /////////////
2128 // Plot
2129 ////////////
2130
02524e30 2131 // Draw histograms for opening angle
2132 TCanvas * copeningangle =new TCanvas("openingangle","Openingangle",800,800);
2133 copeningangle->cd();
70da6c5a 2134 //openingangleppproj->Draw();
2135 //openinganglennproj->Draw("same");
2136 openinganglessproj->Draw();
2137 //openinganglerproj->Draw("same");
02524e30 2138 openingangleosproj->Draw("same");
2139 if(openinganglegammaproj) openinganglegammaproj->Draw("same");
2140 if(openinganglepi0proj) openinganglepi0proj->Draw("same");
70da6c5a 2141 //if(openingangleCproj) openingangleCproj->Draw("same");
2142 //if(openingangleBproj) openingangleBproj->Draw("same");
02524e30 2143 if(openingangleetaproj) openingangleetaproj->Draw("same");
d2af20c5 2144 if(openingangleSplittedTrackssproj) openingangleSplittedTrackssproj->Draw("same");
2145 if(openingangleSplittedTrackosproj) openingangleSplittedTrackosproj->Draw("same");
02524e30 2146 TLegend *lego = new TLegend(0.4,0.6,0.89,0.89);
70da6c5a 2147 //lego->AddEntry(openingangleppproj,"positive-positive","p");
2148 //lego->AddEntry(openinganglennproj,"negative-negative","p");
02524e30 2149 lego->AddEntry(openinganglessproj,"same-sign","p");
70da6c5a 2150 //lego->AddEntry(openinganglerproj,"rotated","p");
02524e30 2151 lego->AddEntry(openingangleosproj,"positive-negative","p");
2152 if(openinganglegammaproj) lego->AddEntry(openinganglegammaproj,"e^{+}e^{-} from #gamma","p");
2153 if(openinganglepi0proj) lego->AddEntry(openinganglepi0proj,"e^{+}e^{-} from #pi^{0}","p");
70da6c5a 2154 //if(openingangleCproj) lego->AddEntry(openingangleCproj,"e^{+}e^{-} from c","p");
2155 //if(openingangleBproj) lego->AddEntry(openingangleBproj,"e^{+}e^{-} from b","p");
02524e30 2156 if(openingangleetaproj) lego->AddEntry(openingangleetaproj,"e^{+}e^{-} from #eta","p");
d2af20c5 2157 if(openingangleSplittedTrackssproj) lego->AddEntry(openingangleSplittedTrackssproj,"Splitted tracks same sign","p");
2158 if(openingangleSplittedTrackosproj) lego->AddEntry(openingangleSplittedTrackosproj,"Splitted tracks opposite sign","p");
02524e30 2159 lego->Draw("same");
70da6c5a 2160
d2af20c5 2161 // Draw histograms for invariant mass
02524e30 2162 TCanvas * cinvmass =new TCanvas("invmass","Invmass",800,800);
2163 cinvmass->cd();
70da6c5a 2164 //invmassppproj->Draw();
2165 //invmassnnproj->Draw("same");
2166 invmassssproj->Draw();
2167 //invmassrproj->Draw("same");
02524e30 2168 invmassosproj->Draw("same");
2169 if(invmassgammaproj) invmassgammaproj->Draw("same");
2170 if(invmasspi0proj) invmasspi0proj->Draw("same");
70da6c5a 2171 //if(invmassCproj) invmassCproj->Draw("same");
2172 //if(invmassBproj) invmassBproj->Draw("same");
02524e30 2173 if(invmassetaproj) invmassetaproj->Draw("same");
d2af20c5 2174 if(invmassSplittedTrackssproj) invmassSplittedTrackssproj->Draw("same");
2175 if(invmassSplittedTrackosproj) invmassSplittedTrackosproj->Draw("same");
02524e30 2176 TLegend *legi = new TLegend(0.4,0.6,0.89,0.89);
70da6c5a 2177 //legi->AddEntry(invmassppproj,"positive-positive","p");
2178 //legi->AddEntry(invmassnnproj,"negative-negative","p");
02524e30 2179 legi->AddEntry(invmassssproj,"same-sign","p");
70da6c5a 2180 //legi->AddEntry(invmassrproj,"rotated","p");
02524e30 2181 legi->AddEntry(invmassosproj,"positive-negative","p");
2182 if(invmassgammaproj) legi->AddEntry(invmassgammaproj,"e^{+}e^{-} from #gamma","p");
2183 if(invmasspi0proj) legi->AddEntry(invmasspi0proj,"e^{+}e^{-} from #pi^{0}","p");
70da6c5a 2184 //if(invmassCproj) legi->AddEntry(invmassCproj,"e^{+}e^{-} from c","p");
2185 //if(invmassBproj) legi->AddEntry(invmassBproj,"e^{+}e^{-} from b","p");
02524e30 2186 if(invmassetaproj) legi->AddEntry(invmassetaproj,"e^{+}e^{-} from #eta","p");
d2af20c5 2187 if(invmassSplittedTrackssproj) legi->AddEntry(invmassSplittedTrackssproj,"Splitted tracks same sign","p");
2188 if(invmassSplittedTrackosproj) legi->AddEntry(invmassSplittedTrackosproj,"Splitted tracks opposite sign","p");
02524e30 2189 legi->Draw("same");
2190
70da6c5a 2191
2192
d2af20c5 2193 // Draw histograms for opening angle 2D
2194 TCanvas * copeningangle2D =new TCanvas("openingangle2D","Openingangle2D",800,800);
2195 copeningangle2D->Divide(6,2);
2196 copeningangle2D->cd(1);
2197 openingangleppproj2D->Draw("lego");
2198 copeningangle2D->cd(2);
2199 openinganglennproj2D->Draw("lego");
2200 copeningangle2D->cd(3);
2201 openinganglessproj2D->Draw("lego");
2202 copeningangle2D->cd(4);
2203 openinganglerproj2D->Draw("lego");
2204 copeningangle2D->cd(5);
2205 openingangleosproj2D->Draw("lego");
2206 copeningangle2D->cd(6);
2207 if(openinganglegammaproj2D) openinganglegammaproj2D->Draw("lego");
2208 copeningangle2D->cd(7);
2209 if(openinganglepi0proj2D) openinganglepi0proj2D->Draw("lego");
2210 copeningangle2D->cd(8);
2211 if(openingangleCproj2D) openingangleCproj2D->Draw("lego");
2212 copeningangle2D->cd(9);
2213 if(openingangleBproj2D) openingangleBproj2D->Draw("lego");
2214 copeningangle2D->cd(10);
2215 if(openingangleetaproj2D) openingangleetaproj2D->Draw("lego");
2216 copeningangle2D->cd(11);
2217 if(openingangleSplittedTrackssproj2D) openingangleSplittedTrackssproj2D->Draw("lego");
2218 copeningangle2D->cd(12);
2219 if(openingangleSplittedTrackosproj2D) openingangleSplittedTrackosproj2D->Draw("lego");
2220
2221 // Draw histograms for invariant mass 2D
2222 TCanvas * cinvmass2D =new TCanvas("invmass2D","Invmass2D",800,800);
2223 cinvmass2D->Divide(6,2);
2224 cinvmass2D->cd(1);
2225 invmassppproj2D->Draw("lego");
2226 cinvmass2D->cd(2);
2227 invmassnnproj2D->Draw("lego");
2228 cinvmass2D->cd(3);
2229 invmassssproj2D->Draw("lego");
2230 cinvmass2D->cd(4);
2231 invmassrproj2D->Draw("lego");
2232 cinvmass2D->cd(5);
2233 invmassosproj2D->Draw("lego");
2234 cinvmass2D->cd(6);
2235 if(invmassgammaproj2D) invmassgammaproj2D->Draw("lego");
2236 cinvmass2D->cd(7);
2237 if(invmasspi0proj2D) invmasspi0proj2D->Draw("lego");
2238 cinvmass2D->cd(8);
2239 if(invmassCproj2D) invmassCproj2D->Draw("lego");
2240 cinvmass2D->cd(9);
2241 if(invmassBproj2D) invmassBproj2D->Draw("lego");
2242 cinvmass2D->cd(10);
2243 if(invmassetaproj2D) invmassetaproj2D->Draw("lego");
2244 cinvmass2D->cd(11);
2245 if(invmassSplittedTrackssproj2D) invmassSplittedTrackssproj2D->Draw("lego");
2246 cinvmass2D->cd(12);
2247 if(invmassSplittedTrackosproj2D) invmassSplittedTrackosproj2D->Draw("lego");
70da6c5a 2248
70da6c5a 2249 /////////////////////////////////////
2250 // Data Radius and chi2Ndf if AliKF
2251 ////////////////////////////////////
2252
2253 TH1F *hDataRadius = dynamic_cast<TH1F *>(fList->FindObject("DataRadius"));
2254 TH1F *hDataChi2Ndf = dynamic_cast<TH1F *>(fList->FindObject("DataChi2Ndf"));
2255
2256 if(hDataRadius || hDataChi2Ndf) {
2257 TCanvas * cDataRadiusChi2Ndf =new TCanvas("CanvasDataRadiusChi2Ndf","CanvasDataRadiusChi2Ndf",800,800);
2258 cDataRadiusChi2Ndf->Divide(2,1);
2259 cDataRadiusChi2Ndf->cd(1);
2260 if(hDataRadius) hDataRadius->Draw();
2261 cDataRadiusChi2Ndf->cd(2);
2262 if(hDataChi2Ndf) hDataChi2Ndf->Draw();
2263 }
2264
2265 ///////////////////////
2266 // Data DCA
2267 //////////////////////
d2af20c5 2268
70da6c5a 2269 TH1F *hDataDCA = dynamic_cast<TH1F *>(fList->FindObject("DataDCA"));
2270
2271 if(hDataDCA) {
2272 TCanvas * cDataDCA =new TCanvas("CanvasDataDCA","CanvasDataDCA",800,800);
2273 cDataDCA->cd(1);
2274 hDataDCA->Draw();
2275 }
2276
2277 /////////////////////////////////////
2278 // MC Radius and chi2Ndf if AliKF
2279 ////////////////////////////////////
2280
2281 TH2F *hMCRadius = dynamic_cast<TH2F *>(fList->FindObject("MCRadius"));
2282 TH2F *hMCChi2Ndf = dynamic_cast<TH2F *>(fList->FindObject("MCChi2Ndf"));
2283
2284 if(hMCRadius || hMCChi2Ndf) {
2285 TCanvas * cMCRadiusChi2Ndf =new TCanvas("CanvasMCRadiusChi2Ndf","CanvasMCRadiusChi2Ndf",800,800);
2286 cMCRadiusChi2Ndf->Divide(2,1);
2287 cMCRadiusChi2Ndf->cd(1);
2288 //TH1D *hMCRadiusBackground = hMCRadius->ProjectionX("MCRadiusBackGround",1,1,"e");
2289 TH1D *hMCRadiusGamma = hMCRadius->ProjectionX("MCRadiusGamma",2,2,"e");
2290 TH1D *hMCRadiusPi0 = hMCRadius->ProjectionX("MCRadiusPi0",3,3,"e");
2291 TH1D *hMCRadiusEta = hMCRadius->ProjectionX("MCRadiusEta",4,4,"e");
2292 TH1D *hMCRadiusC = hMCRadius->ProjectionX("MCRadiusC",5,5,"e");
2293 TH1D *hMCRadiusB = hMCRadius->ProjectionX("MCRadiusB",6,6,"e");
2294 //hMCRadiusBackground->Draw();
2295 hMCRadiusGamma->Draw();
2296 hMCRadiusPi0->Draw("same");
2297 hMCRadiusEta->Draw("same");
2298 hMCRadiusC->Draw("same");
2299 hMCRadiusB->Draw("same");
2300 TLegend *legRadius = new TLegend(0.4,0.6,0.89,0.89);
2301 //legRadius->AddEntry(hMCRadiusBackground,"Background","p");
2302 legRadius->AddEntry(hMCRadiusGamma,"#gamma","p");
2303 legRadius->AddEntry(hMCRadiusPi0,"#pi^{0}","p");
2304 legRadius->AddEntry(hMCRadiusEta,"#eta","p");
2305 legRadius->AddEntry(hMCRadiusC,"c","p");
2306 legRadius->AddEntry(hMCRadiusB,"b","p");
2307 legRadius->Draw("same");
2308 cMCRadiusChi2Ndf->cd(2);
2309 //TH1D *hMCChi2NdfBackground = hMCChi2Ndf->ProjectionX("MCChi2NdfBackGround",1,1,"e");
2310 TH1D *hMCChi2NdfGamma = hMCChi2Ndf->ProjectionX("MCChi2NdfGamma",2,2,"e");
2311 TH1D *hMCChi2NdfPi0 = hMCChi2Ndf->ProjectionX("MCChi2NdfPi0",3,3,"e");
2312 TH1D *hMCChi2NdfEta = hMCChi2Ndf->ProjectionX("MCChi2NdfEta",4,4,"e");
2313 TH1D *hMCChi2NdfC = hMCChi2Ndf->ProjectionX("MCChi2NdfC",5,5,"e");
2314 TH1D *hMCChi2NdfB = hMCChi2Ndf->ProjectionX("MCChi2NdfB",6,6,"e");
2315 //hMCChi2NdfBackground->Draw();
2316 hMCChi2NdfGamma->Draw();
2317 hMCChi2NdfPi0->Draw("same");
2318 hMCChi2NdfEta->Draw("same");
2319 hMCChi2NdfC->Draw("same");
2320 hMCChi2NdfB->Draw("same");
2321 TLegend *legChi2Ndf = new TLegend(0.4,0.6,0.89,0.89);
2322 //legChi2Ndf->AddEntry(hMCChi2NdfBackground,"Background","p");
2323 legChi2Ndf->AddEntry(hMCChi2NdfGamma,"#gamma","p");
2324 legChi2Ndf->AddEntry(hMCChi2NdfPi0,"#pi^{0}","p");
2325 legChi2Ndf->AddEntry(hMCChi2NdfEta,"#eta","p");
2326 legChi2Ndf->AddEntry(hMCChi2NdfC,"c","p");
2327 legChi2Ndf->AddEntry(hMCChi2NdfB,"b","p");
2328 legChi2Ndf->Draw("same");
2329 }
2330
2331 ///////////////////////
2332 // MC DCA
2333 //////////////////////
2334
2335 TH2F *hMCDCA = dynamic_cast<TH2F *>(fList->FindObject("MCDCA"));
2336
2337 if(hMCDCA) {
2338 TCanvas * cMCDCA =new TCanvas("CanvasMCDCA","CanvasMCDCA",800,800);
2339 cMCDCA->cd(1);
2340 //TH1D *hMCDCABackground = hMCDCA->ProjectionX("MCDCABackGround",1,1,"e");
2341 TH1D *hMCDCAGamma = hMCDCA->ProjectionX("MCDCAGamma",2,2,"e");
2342 TH1D *hMCDCAPi0 = hMCDCA->ProjectionX("MCDCAPi0",3,3,"e");
2343 TH1D *hMCDCAEta = hMCDCA->ProjectionX("MCDCAEta",4,4,"e");
2344 TH1D *hMCDCAC = hMCDCA->ProjectionX("MCDCAC",5,5,"e");
2345 TH1D *hMCDCAB = hMCDCA->ProjectionX("MCDCAB",6,6,"e");
2346 //hMCDCABackground->Draw();
2347 hMCDCAGamma->Draw();
2348 hMCDCAPi0->Draw("same");
2349 hMCDCAEta->Draw("same");
2350 hMCDCAC->Draw("same");
2351 hMCDCAB->Draw("same");
2352 TLegend *legDCA = new TLegend(0.4,0.6,0.89,0.89);
2353 //legDCA->AddEntry(hMCDCABackground,"Background","p");
2354 legDCA->AddEntry(hMCDCAGamma,"#gamma","p");
2355 legDCA->AddEntry(hMCDCAPi0,"#pi^{0}","p");
2356 legDCA->AddEntry(hMCDCAEta,"#eta","p");
2357 legDCA->AddEntry(hMCDCAC,"c","p");
2358 legDCA->AddEntry(hMCDCAB,"b","p");
2359 legDCA->Draw("same");
2360 }
2361
2362
2363 /////////////////////////
2364 // PID Partner Signal
2365 /////////////////////////
2366 TF1 *betheBlochElectron = 0x0;
2367 TF1 *betheBlochMuon = 0x0;
2368 TF1 *betheBlochPion = 0x0;
2369 TF1 *betheBlochKaon = 0x0;
2370 TF1 *betheBlochProton = 0x0;
2371
2372 TH2F *hsignalPidPartner0 = dynamic_cast<TH2F *>(fList->FindObject("TPCPartner0"));
2373 TH2F *hsignalPidPartner1 = dynamic_cast<TH2F *>(fList->FindObject("TPCPartner1"));
2374 if((!hsignalPidPartner0) && (!hsignalPidPartner1)) {
2375 hsignalPidPartner0 = dynamic_cast<TH2F *>(fList->FindObject("ITSPartner0"));
2376 hsignalPidPartner1 = dynamic_cast<TH2F *>(fList->FindObject("ITSPartner1"));
2377
2378 betheBlochElectron = new TF1("betheBlochElectron",BetheBlochElectronITS,0.1,10.0,0);
2379 betheBlochMuon = new TF1("betheBlochMuon",BetheBlochMuonITS,0.1,10.0,0);
2380 betheBlochPion = new TF1("betheBlochPion",BetheBlochPionITS,0.1,10.0,0);
2381 betheBlochKaon = new TF1("betheBlochKaon",BetheBlochKaonITS,0.1,10.0,0);
2382 betheBlochProton = new TF1("betheBlochProton",BetheBlochProtonITS,0.1,10.0,0);
2383
2384 }
2385 else {
2386
2387 betheBlochElectron = new TF1("betheBlochElectron",BetheBlochElectronTPC,0.1,10.0,0);
2388 betheBlochMuon = new TF1("betheBlochMuon",BetheBlochMuonTPC,0.1,10.0,0);
2389 betheBlochPion = new TF1("betheBlochPion",BetheBlochPionTPC,0.1,10.0,0);
2390 betheBlochKaon = new TF1("betheBlochKaon",BetheBlochKaonTPC,0.1,10.0,0);
2391 betheBlochProton = new TF1("betheBlochProton",BetheBlochProtonTPC,0.1,10.0,0);
2392
2393 }
2394
2395
2396 if((hsignalPidPartner0) || (hsignalPidPartner1)) {
2397 TCanvas * cPidSignal =new TCanvas("cPidSignal","cPidSignal",800,800);
2398 cPidSignal->Divide(2,1);
2399 cPidSignal->cd(1);
2400 if(hsignalPidPartner0) hsignalPidPartner0->Draw("colz");
2401 if(betheBlochElectron) betheBlochElectron->Draw("same");
2402 if(betheBlochMuon) betheBlochMuon->Draw("same");
2403 if(betheBlochPion) betheBlochPion->Draw("same");
2404 if(betheBlochKaon) betheBlochKaon->Draw("same");
2405 if(betheBlochProton) betheBlochProton->Draw("same");
2406 cPidSignal->cd(2);
2407 if(hsignalPidPartner1) hsignalPidPartner1->Draw("colz");
2408 if(betheBlochElectron) betheBlochElectron->Draw("same");
2409 if(betheBlochMuon) betheBlochMuon->Draw("same");
2410 if(betheBlochPion) betheBlochPion->Draw("same");
2411 if(betheBlochKaon) betheBlochKaon->Draw("same");
2412 if(betheBlochProton) betheBlochProton->Draw("same");
2413 }
2414
2415 THnSparseF *hsSparseITSsignal = dynamic_cast<THnSparseF *>(fList->FindObject("SparseITSsignal"));
2416 if(hsSparseITSsignal) {
2417
2418
2419 TH2D *sddsdd = hsSparseITSsignal->Projection(1,2);
2420 TH2D *ssdssd = hsSparseITSsignal->Projection(3,4);
2421 TH2D *sddssda = hsSparseITSsignal->Projection(1,3);
2422 TH2D *sddssdb = hsSparseITSsignal->Projection(2,4);
2423 TH2D *sddssdc = hsSparseITSsignal->Projection(1,4);
2424 TH2D *sddssdd = hsSparseITSsignal->Projection(2,3);
2425
2426 TCanvas * cITSSignal =new TCanvas("cITSSignal","cITSSignal",800,800);
2427 cITSSignal->Divide(2,3);
2428 cITSSignal->cd(1);
2429 sddsdd->Draw("colz");
2430 cITSSignal->cd(2);
2431 ssdssd->Draw("colz");
2432 cITSSignal->cd(3);
2433 sddssda->Draw("colz");
2434 cITSSignal->cd(4);
2435 sddssdb->Draw("colz");
2436 cITSSignal->cd(5);
2437 sddssdc->Draw("colz");
2438 cITSSignal->cd(6);
2439 sddssdd->Draw("colz");
2440
2441 }
2442
2443 THnSparseF *hsSparseITSsignalSplit = dynamic_cast<THnSparseF *>(fList->FindObject("SparseITSsignalSplit"));
2444 if(hsSparseITSsignalSplit) {
2445
2446 // no splitted
2447 hsSparseITSsignalSplit->GetAxis(0)->SetRange(1,1);
2448
2449 TH1D *layerITS2 = hsSparseITSsignalSplit->Projection(1);
2450 TH1D *layerITS3 = hsSparseITSsignalSplit->Projection(2);
2451 TH1D *layerITS4 = hsSparseITSsignalSplit->Projection(3);
2452 TH1D *layerITS5 = hsSparseITSsignalSplit->Projection(4);
2453
2454 // splitted
2455 hsSparseITSsignalSplit->GetAxis(0)->SetRange(2,2);
2456
2457 TH1D *layerITS2s = hsSparseITSsignalSplit->Projection(1);
2458 TH1D *layerITS3s = hsSparseITSsignalSplit->Projection(2);
2459 TH1D *layerITS4s = hsSparseITSsignalSplit->Projection(3);
2460 TH1D *layerITS5s = hsSparseITSsignalSplit->Projection(4);
2461
2462 TCanvas * cITSSignalSplit =new TCanvas("cITSSignalSplit","cITSSignalSplit",800,800);
2463 cITSSignalSplit->Divide(2,2);
2464 cITSSignalSplit->cd(1);
2465 layerITS2->Draw();
2466 layerITS2s->Draw("same");
2467 TLegend *legITS2 = new TLegend(0.4,0.6,0.89,0.89);
2468 legITS2->AddEntry(layerITS2,"No splitted","p");
2469 legITS2->AddEntry(layerITS2s,"Splitted","p");
2470 legITS2->Draw("same");
2471 cITSSignalSplit->cd(2);
2472 layerITS3->Draw();
2473 layerITS3s->Draw("same");
2474 TLegend *legITS3 = new TLegend(0.4,0.6,0.89,0.89);
2475 legITS3->AddEntry(layerITS3,"No splitted","p");
2476 legITS3->AddEntry(layerITS3s,"Splitted","p");
2477 legITS3->Draw("same");
2478 cITSSignalSplit->cd(3);
2479 layerITS4->Draw();
2480 layerITS4s->Draw("same");
2481 TLegend *legITS4 = new TLegend(0.4,0.6,0.89,0.89);
2482 legITS4->AddEntry(layerITS4,"No splitted","p");
2483 legITS4->AddEntry(layerITS4s,"Splitted","p");
2484 legITS4->Draw("same");
2485 cITSSignalSplit->cd(4);
2486 layerITS5->Draw();
2487 layerITS5s->Draw("same");
2488 TLegend *legITS5 = new TLegend(0.4,0.6,0.89,0.89);
2489 legITS5->AddEntry(layerITS5,"No splitted","p");
2490 legITS5->AddEntry(layerITS5s,"Splitted","p");
2491 legITS5->Draw("same");
2492
2493
e3ae862b 2494 }
70da6c5a 2495
e3ae862b 2496 ////////////////////////
2497 // Cut efficiencies
2498 ////////////////////////
2499
2500 THnSparseF *hsSparseMCe = dynamic_cast<THnSparseF *>(fList->FindObject("CutPassedMC"));
2501 if(!hsSparseMCe) return;
2502
2503 // init histos
2504 TAxis *axissources = hsSparseMCe->GetAxis(2);
2505 Int_t nbsources = axissources->GetNbins();
2506 TAxis *axiscuts = hsSparseMCe->GetAxis(1);
2507 Int_t nbcuts = axiscuts->GetNbins();
2508 TH1D **histopassedcuts = new TH1D*[nbsources*nbcuts];
2509 Double_t *nbEntriesCuts = new Double_t[nbsources*nbcuts];
2510 for(Int_t k =0; k < nbsources*nbcuts; k++){
2511 nbEntriesCuts[k] = 0.0;
2512 histopassedcuts[k] = 0x0;
2513 }
2514
2515 //printf("Number of cuts %d\n",nbcuts);
2516
2517 // canvas
2518 TCanvas * chsSparseMCeeff =new TCanvas("hsSparseMCeeffDebug","hsSparseMCeeffDebug",800,800);
2519 chsSparseMCeeff->Divide(3,1);
2520
2521 // histos
2522 for(Int_t sourceid = 0; sourceid < nbsources; sourceid++) {
2523 hsSparseMCe->GetAxis(2)->SetRange(sourceid+1,sourceid+1);
2524 for(Int_t cut = 0; cut < nbcuts; cut++){
2525 hsSparseMCe->GetAxis(1)->SetRange(cut+1,cut+1);
2526 histopassedcuts[sourceid*nbcuts+cut] = hsSparseMCe->Projection(0);
2527 hsSparseMCe->GetAxis(1)->SetRange(1,hsSparseMCe->GetAxis(1)->GetNbins());
2528 }
2529 hsSparseMCe->GetAxis(2)->SetRange(1,hsSparseMCe->GetAxis(2)->GetNbins());
2530 }
2531
2532 // calcul efficiencies
2533
2534 // histos
2535 for(Int_t sourceid = 0; sourceid < nbsources; sourceid++) {
2536 // Next is compared to the partner tracked
2537 for(Int_t cut = 2; cut < nbcuts; cut++){
2538 nbEntriesCuts[sourceid*nbcuts+cut] = histopassedcuts[sourceid*nbcuts+cut]->GetEntries();
2539 if(histopassedcuts[sourceid*nbcuts+1]->GetEntries() > 0.0) histopassedcuts[sourceid*nbcuts+cut]->Divide(histopassedcuts[sourceid*nbcuts+1]);
2540 }
2541 // First one is if the partner is tracked.
2542 nbEntriesCuts[sourceid*nbcuts+1] = histopassedcuts[sourceid*nbcuts+1]->GetEntries();
2543 if(histopassedcuts[sourceid*nbcuts]->GetEntries() > 0.0) histopassedcuts[sourceid*nbcuts+1]->Divide(histopassedcuts[sourceid*nbcuts]);
2544 // First one is input
2545 nbEntriesCuts[sourceid*nbcuts] = histopassedcuts[sourceid*nbcuts]->GetEntries();
2546 }
2547
2548 // ratios
2549 for(Int_t sourceid = 0; sourceid < nbsources; sourceid++) {
2550 for(Int_t cut = 1; cut < nbcuts; cut++){
2551 if(nbEntriesCuts[sourceid*nbcuts] > 0.0) nbEntriesCuts[sourceid*nbcuts+cut] = nbEntriesCuts[sourceid*nbcuts+cut]/nbEntriesCuts[sourceid*nbcuts];
2552 }
2553 }
2554 TH1F *ratioHistoEntriesGamma = new TH1F("ratioHistoEntriesGamma","", nbcuts-1, 0.0, nbcuts-1.0);
2555 TH1F *ratioHistoEntriesPi0 = new TH1F("ratioHistoEntriesPi0","", nbcuts-1, 0.0, nbcuts-1.0);
2556 TH1F *ratioHistoEntriesC = new TH1F("ratioHistoEntriesC","", nbcuts-1, 0.0, nbcuts-1.0);
2557 for(Int_t k = 1; k < nbcuts; k++){
2558 if((nbcuts+k) < (nbsources*nbcuts)) ratioHistoEntriesGamma->SetBinContent(k,nbEntriesCuts[nbcuts+k]);
2559 if((2*nbcuts+k) < (nbsources*nbcuts)) ratioHistoEntriesPi0->SetBinContent(k,nbEntriesCuts[2*nbcuts+k]);
2560 if((4*nbcuts+k) < (nbsources*nbcuts)) ratioHistoEntriesC->SetBinContent(k,nbEntriesCuts[4*nbcuts+k]);
2561 }
2562 //
2563 TAxis *xAxisGamma = ratioHistoEntriesGamma->GetXaxis();
2564 xAxisGamma->SetBinLabel(1,"Partner tracked");
2565 xAxisGamma->SetBinLabel(2,"Opposite sign");
2566 xAxisGamma->SetBinLabel(3,"Single Track Cut");
2567 xAxisGamma->SetBinLabel(4,"Shared Clusters");
2568 xAxisGamma->SetBinLabel(5,"PID");
2569 xAxisGamma->SetBinLabel(6,"DCA");
2570 xAxisGamma->SetBinLabel(7,"Chi^{2}/Ndf");
2571 xAxisGamma->SetBinLabel(8,"Opening angle");
2572 xAxisGamma->SetBinLabel(9,"Invariant mass");
2573 //
2574 TAxis *xAxisPi0 = ratioHistoEntriesPi0->GetXaxis();
2575 xAxisPi0->SetBinLabel(1,"Partner tracked");
2576 xAxisPi0->SetBinLabel(2,"Opposite sign");
2577 xAxisPi0->SetBinLabel(3,"Single Track Cut");
2578 xAxisPi0->SetBinLabel(4,"Shared Clusters");
2579 xAxisPi0->SetBinLabel(5,"PID");
2580 xAxisPi0->SetBinLabel(6,"DCA");
2581 xAxisPi0->SetBinLabel(7,"Chi^{2}/Ndf");
2582 xAxisPi0->SetBinLabel(8,"Opening angle");
2583 xAxisPi0->SetBinLabel(9,"Invariant mass");
2584 //
2585 TAxis *xAxisC = ratioHistoEntriesC->GetXaxis();
2586 xAxisC->SetBinLabel(1,"Partner tracked");
2587 xAxisC->SetBinLabel(2,"Opposite sign");
2588 xAxisC->SetBinLabel(3,"Single Track Cut");
2589 xAxisC->SetBinLabel(4,"Shared Clusters");
2590 xAxisC->SetBinLabel(5,"PID");
2591 xAxisC->SetBinLabel(6,"DCA");
2592 xAxisC->SetBinLabel(7,"Chi^{2}/Ndf");
2593 xAxisC->SetBinLabel(8,"Opening angle");
2594 xAxisC->SetBinLabel(9,"Invariant mass");
2595 //
2596 TCanvas * cRatioHistoEntries =new TCanvas("cRatioHistoEntries","cRatioHistoEntries",800,800);
2597 cRatioHistoEntries->cd(1);
2598 ratioHistoEntriesGamma->SetStats(0);
2599 ratioHistoEntriesGamma->Draw();
2600 ratioHistoEntriesPi0->SetStats(0);
2601 ratioHistoEntriesPi0->Draw("same");
2602 ratioHistoEntriesC->SetStats(0);
2603 //ratioHistoEntriesC->Draw("same");
2604 TLegend *legEntries = new TLegend(0.4,0.6,0.89,0.89);
2605 legEntries->AddEntry(ratioHistoEntriesGamma,"#gamma","l");
2606 legEntries->AddEntry(ratioHistoEntriesPi0,"#pi^{0}","l");
2607 //legEntries->AddEntry(ratioHistoEntriesC,"c","p");
2608 legEntries->Draw("same");
2609
2610 // plot Debug
2611 Int_t source = 1;
2612 chsSparseMCeeff->cd(1);
ccc37cdc 2613 if(((source*nbcuts+0)> (nbsources*nbcuts-1)) || ((source*nbcuts+1)> (nbsources*nbcuts-1)) || ((source*nbcuts+2)> (nbsources*nbcuts-1)) || ((source*nbcuts+3)> (nbsources*nbcuts-1)) || ((source*nbcuts+4)> (nbsources*nbcuts-1)) || ((source*nbcuts+5)> (nbsources*nbcuts-1)) || ((source*nbcuts+6)> (nbsources*nbcuts-1)) || ((source*nbcuts+7)> (nbsources*nbcuts-1)) || ((source*nbcuts+8)> (nbsources*nbcuts-1)) || ((source*nbcuts+9)> (nbsources*nbcuts-1))) {
2614 delete [] histopassedcuts;
2615 delete [] nbEntriesCuts;
2616 return;
2617 }
2618 if((!histopassedcuts[source*nbcuts+0]) || (!histopassedcuts[source*nbcuts+1]) || (!histopassedcuts[source*nbcuts+2]) || (!histopassedcuts[source*nbcuts+3]) || (!histopassedcuts[source*nbcuts+4]) || (!histopassedcuts[source*nbcuts+5]) || (!histopassedcuts[source*nbcuts+6]) || (!histopassedcuts[source*nbcuts+7]) || (!histopassedcuts[source*nbcuts+8]) || (!histopassedcuts[source*nbcuts+9])) {
2619 delete [] histopassedcuts;
2620 delete [] nbEntriesCuts;
2621 return;
2622 }
e3ae862b 2623 histopassedcuts[source*nbcuts+0]->SetTitle("#gamma");
2624 histopassedcuts[source*nbcuts+1]->SetTitle("#gamma");
2625 histopassedcuts[source*nbcuts+2]->SetTitle("#gamma");
2626 histopassedcuts[source*nbcuts+3]->SetTitle("#gamma");
2627 histopassedcuts[source*nbcuts+4]->SetTitle("#gamma");
2628 histopassedcuts[source*nbcuts+5]->SetTitle("#gamma");
2629 histopassedcuts[source*nbcuts+6]->SetTitle("#gamma");
2630 histopassedcuts[source*nbcuts+7]->SetTitle("#gamma");
2631 histopassedcuts[source*nbcuts+8]->SetTitle("#gamma");
2632 histopassedcuts[source*nbcuts+9]->SetTitle("#gamma");
2633 //histopassedcuts[source*nbcuts+0]->SetStats(0);
2634 histopassedcuts[source*nbcuts+1]->SetStats(0);
2635 histopassedcuts[source*nbcuts+2]->SetStats(0);
2636 histopassedcuts[source*nbcuts+3]->SetStats(0);
2637 histopassedcuts[source*nbcuts+4]->SetStats(0);
2638 histopassedcuts[source*nbcuts+5]->SetStats(0);
2639 histopassedcuts[source*nbcuts+6]->SetStats(0);
2640 histopassedcuts[source*nbcuts+7]->SetStats(0);
2641 histopassedcuts[source*nbcuts+8]->SetStats(0);
2642 histopassedcuts[source*nbcuts+9]->SetStats(0);
2643 //histopassedcuts[source*nbcuts+0]->Draw();
2644 //histopassedcuts[source*nbcuts+1]->Draw("");
2645 histopassedcuts[source*nbcuts+2]->Draw();
2646 histopassedcuts[source*nbcuts+3]->Draw("same");
2647 //histopassedcuts[source*nbcuts+4]->Draw("same");
2648 histopassedcuts[source*nbcuts+5]->Draw("same");
2649 histopassedcuts[source*nbcuts+6]->Draw("same");
2650 //histopassedcuts[source*nbcuts+7]->Draw("same");
2651 histopassedcuts[source*nbcuts+8]->Draw("same");
2652 histopassedcuts[source*nbcuts+9]->Draw("same");
2653 TLegend *legb = new TLegend(0.4,0.6,0.89,0.89);
2654 //legb->AddEntry(histopassedcuts[source*nbcuts+0],"all","p");
2655 //legb->AddEntry(histopassedcuts[source*nbcuts+1],"Partner tracked","p");
2656 legb->AddEntry(histopassedcuts[source*nbcuts+2],"Opposite sign","p");
2657 legb->AddEntry(histopassedcuts[source*nbcuts+3],"SingleTrackPart","p");
2658 //legb->AddEntry(histopassedcuts[source*nbcuts+4],"SharedCluster","p");
2659 legb->AddEntry(histopassedcuts[source*nbcuts+5],"PID","p");
2660 legb->AddEntry(histopassedcuts[source*nbcuts+6],"DCA","p");
2661 //legb->AddEntry(histopassedcuts[source*nbcuts+7],"Chi2Ndf","p");
2662 legb->AddEntry(histopassedcuts[source*nbcuts+8],"OpeningAngle","p");
2663 legb->AddEntry(histopassedcuts[source*nbcuts+9],"InvMass","p");
2664 legb->Draw("same");
2665
2666 source = 2;
ccc37cdc 2667 if(((source*nbcuts+0)> (nbsources*nbcuts-1)) || ((source*nbcuts+1)> (nbsources*nbcuts-1)) || ((source*nbcuts+2)> (nbsources*nbcuts-1)) || ((source*nbcuts+3)> (nbsources*nbcuts-1)) || ((source*nbcuts+4)> (nbsources*nbcuts-1)) || ((source*nbcuts+5)> (nbsources*nbcuts-1)) || ((source*nbcuts+6)> (nbsources*nbcuts-1)) || ((source*nbcuts+7)> (nbsources*nbcuts-1)) || ((source*nbcuts+8)> (nbsources*nbcuts-1)) || ((source*nbcuts+9)> (nbsources*nbcuts-1))) {
2668 delete [] histopassedcuts;
2669 delete [] nbEntriesCuts;
2670 return;
2671 }
2672 if((!histopassedcuts[source*nbcuts+0]) || (!histopassedcuts[source*nbcuts+1]) || (!histopassedcuts[source*nbcuts+2]) || (!histopassedcuts[source*nbcuts+3]) || (!histopassedcuts[source*nbcuts+4]) || (!histopassedcuts[source*nbcuts+5]) || (!histopassedcuts[source*nbcuts+6]) || (!histopassedcuts[source*nbcuts+7]) || (!histopassedcuts[source*nbcuts+8]) || (!histopassedcuts[source*nbcuts+9])) {
2673 delete [] histopassedcuts;
2674 delete [] nbEntriesCuts;
2675 return;
2676 }
e3ae862b 2677 chsSparseMCeeff->cd(2);
2678 histopassedcuts[source*nbcuts+0]->SetTitle("#pi^{0}");
2679 histopassedcuts[source*nbcuts+1]->SetTitle("#pi^{0}");
2680 histopassedcuts[source*nbcuts+2]->SetTitle("#pi^{0}");
2681 histopassedcuts[source*nbcuts+3]->SetTitle("#pi^{0}");
2682 histopassedcuts[source*nbcuts+4]->SetTitle("#pi^{0}");
2683 histopassedcuts[source*nbcuts+5]->SetTitle("#pi^{0}");
2684 histopassedcuts[source*nbcuts+6]->SetTitle("#pi^{0}");
2685 histopassedcuts[source*nbcuts+7]->SetTitle("#pi^{0}");
2686 histopassedcuts[source*nbcuts+8]->SetTitle("#pi^{0}");
2687 histopassedcuts[source*nbcuts+9]->SetTitle("#pi^{0}");
2688 //histopassedcuts[source*nbcuts+0]->SetStats(0);
2689 histopassedcuts[source*nbcuts+1]->SetStats(0);
2690 histopassedcuts[source*nbcuts+2]->SetStats(0);
2691 histopassedcuts[source*nbcuts+3]->SetStats(0);
2692 histopassedcuts[source*nbcuts+4]->SetStats(0);
2693 histopassedcuts[source*nbcuts+5]->SetStats(0);
2694 histopassedcuts[source*nbcuts+6]->SetStats(0);
2695 histopassedcuts[source*nbcuts+7]->SetStats(0);
2696 histopassedcuts[source*nbcuts+8]->SetStats(0);
2697 histopassedcuts[source*nbcuts+9]->SetStats(0);
2698 //histopassedcuts[source*nbcuts+0]->Draw();
2699 //histopassedcuts[source*nbcuts+1]->Draw();
2700 histopassedcuts[source*nbcuts+2]->Draw();
2701 histopassedcuts[source*nbcuts+3]->Draw("same");
2702 //histopassedcuts[source*nbcuts+4]->Draw("same");
2703 histopassedcuts[source*nbcuts+5]->Draw("same");
2704 histopassedcuts[source*nbcuts+6]->Draw("same");
2705 //histopassedcuts[source*nbcuts+7]->Draw("same");
2706 histopassedcuts[source*nbcuts+8]->Draw("same");
2707 histopassedcuts[source*nbcuts+9]->Draw("same");
2708 TLegend *legc = new TLegend(0.4,0.6,0.89,0.89);
2709 //legc->AddEntry(histopassedcuts[source*nbcuts+0],"all","p");
2710 //legc->AddEntry(histopassedcuts[source*nbcuts+1],"Partner tracked","p");
2711 legc->AddEntry(histopassedcuts[source*nbcuts+2],"Opposite sign","p");
2712 legc->AddEntry(histopassedcuts[source*nbcuts+3],"SingleTrackPart","p");
2713 //legc->AddEntry(histopassedcuts[source*nbcuts+4],"SharedCluster","p");
2714 legc->AddEntry(histopassedcuts[source*nbcuts+5],"PID","p");
2715 legc->AddEntry(histopassedcuts[source*nbcuts+6],"DCA","p");
2716 //legc->AddEntry(histopassedcuts[source*nbcuts+7],"Chi2Ndf","p");
2717 legc->AddEntry(histopassedcuts[source*nbcuts+8],"OpeningAngle","p");
2718 legc->AddEntry(histopassedcuts[source*nbcuts+9],"InvMass","p");
2719 legc->Draw("same");
2720
2721 source = 4;
ccc37cdc 2722 if(((source*nbcuts+0)> (nbsources*nbcuts-1)) || ((source*nbcuts+1)> (nbsources*nbcuts-1)) || ((source*nbcuts+2)> (nbsources*nbcuts-1)) || ((source*nbcuts+3)> (nbsources*nbcuts-1)) || ((source*nbcuts+4)> (nbsources*nbcuts-1)) || ((source*nbcuts+5)> (nbsources*nbcuts-1)) || ((source*nbcuts+6)> (nbsources*nbcuts-1)) || ((source*nbcuts+7)> (nbsources*nbcuts-1)) || ((source*nbcuts+8)> (nbsources*nbcuts-1)) || ((source*nbcuts+9)> (nbsources*nbcuts-1))) {
2723 delete [] histopassedcuts;
2724 delete [] nbEntriesCuts;
2725 return;
2726 }
2727 if((!histopassedcuts[source*nbcuts+0]) || (!histopassedcuts[source*nbcuts+1]) || (!histopassedcuts[source*nbcuts+2]) || (!histopassedcuts[source*nbcuts+3]) || (!histopassedcuts[source*nbcuts+4]) || (!histopassedcuts[source*nbcuts+5]) || (!histopassedcuts[source*nbcuts+6]) || (!histopassedcuts[source*nbcuts+7]) || (!histopassedcuts[source*nbcuts+8]) || (!histopassedcuts[source*nbcuts+9])) {
2728 delete [] histopassedcuts;
2729 delete [] nbEntriesCuts;
2730 return;
2731 }
e3ae862b 2732 chsSparseMCeeff->cd(3);
2733 histopassedcuts[source*nbcuts+0]->SetTitle("C");
2734 histopassedcuts[source*nbcuts+1]->SetTitle("C");
2735 histopassedcuts[source*nbcuts+2]->SetTitle("C");
2736 histopassedcuts[source*nbcuts+3]->SetTitle("C");
2737 histopassedcuts[source*nbcuts+4]->SetTitle("C");
2738 histopassedcuts[source*nbcuts+5]->SetTitle("C");
2739 histopassedcuts[source*nbcuts+6]->SetTitle("C");
2740 histopassedcuts[source*nbcuts+7]->SetTitle("C");
2741 histopassedcuts[source*nbcuts+8]->SetTitle("C");
2742 histopassedcuts[source*nbcuts+9]->SetTitle("C");
2743 //histopassedcuts[source*nbcuts+0]->SetStats(0);
2744 histopassedcuts[source*nbcuts+1]->SetStats(0);
2745 histopassedcuts[source*nbcuts+2]->SetStats(0);
2746 histopassedcuts[source*nbcuts+3]->SetStats(0);
2747 histopassedcuts[source*nbcuts+4]->SetStats(0);
2748 histopassedcuts[source*nbcuts+5]->SetStats(0);
2749 histopassedcuts[source*nbcuts+6]->SetStats(0);
2750 histopassedcuts[source*nbcuts+7]->SetStats(0);
2751 histopassedcuts[source*nbcuts+8]->SetStats(0);
2752 histopassedcuts[source*nbcuts+9]->SetStats(0);
2753 //histopassedcuts[source*nbcuts+0]->Draw();
2754 //histopassedcuts[source*nbcuts+1]->Draw();
2755 histopassedcuts[source*nbcuts+2]->Draw();
2756 histopassedcuts[source*nbcuts+3]->Draw("same");
2757 //histopassedcuts[source*nbcuts+4]->Draw("same");
2758 histopassedcuts[source*nbcuts+5]->Draw("same");
2759 histopassedcuts[source*nbcuts+6]->Draw("same");
2760 //histopassedcuts[source*nbcuts+7]->Draw("same");
2761 histopassedcuts[source*nbcuts+8]->Draw("same");
2762 histopassedcuts[source*nbcuts+9]->Draw("same");
2763 TLegend *lege = new TLegend(0.4,0.6,0.89,0.89);
2764 //lege->AddEntry(histopassedcuts[source*nbcuts+0],"all","p");
2765 //lege->AddEntry(histopassedcuts[source*nbcuts+1],"Partner tracked","p");
2766 lege->AddEntry(histopassedcuts[source*nbcuts+2],"Opposite sign","p");
2767 lege->AddEntry(histopassedcuts[source*nbcuts+3],"SingleTrackPart","p");
2768 //lege->AddEntry(histopassedcuts[source*nbcuts+4],"SharedCluster","p");
2769 lege->AddEntry(histopassedcuts[source*nbcuts+5],"PID","p");
2770 lege->AddEntry(histopassedcuts[source*nbcuts+6],"DCA","p");
2771 //lege->AddEntry(histopassedcuts[source*nbcuts+7],"Chi2Ndf","p");
2772 lege->AddEntry(histopassedcuts[source*nbcuts+8],"OpeningAngle","p");
2773 lege->AddEntry(histopassedcuts[source*nbcuts+9],"InvMass","p");
2774 lege->Draw("same");
2775
2776 //////////////////////
2777 // Input
2778 //////////////////////
2779
2780 TCanvas * chsSparseMCein =new TCanvas("hsSparseMCeinput","hsSparseMCeinput",800,800);
2781 chsSparseMCein->cd(1);
2782 Double_t nbGamma = 0.0;
2783 source = 1;
ccc37cdc 2784 if(((source*nbcuts+0)> (nbsources*nbcuts-1)) || (!histopassedcuts[source*nbcuts+0])) {
2785 delete [] histopassedcuts;
2786 delete [] nbEntriesCuts;
2787 return;
2788 }
e3ae862b 2789 nbGamma = histopassedcuts[source*nbcuts+0]->GetEntries();
2790 histopassedcuts[source*nbcuts+0]->SetStats(0);
2791 histopassedcuts[source*nbcuts+0]->Draw();
2792 TLegend *leginput = new TLegend(0.4,0.6,0.89,0.89);
2793 leginput->AddEntry(histopassedcuts[source*nbcuts+0],"#gamma","p");
2794 Double_t nbPi0 = 0.0;
2795 source = 2;
ccc37cdc 2796 if(((source*nbcuts+0)> (nbsources*nbcuts-1)) || (!histopassedcuts[source*nbcuts+0])) {
2797 delete [] histopassedcuts;
2798 delete [] nbEntriesCuts;
2799 return;
2800 }
e3ae862b 2801 nbPi0 = histopassedcuts[source*nbcuts+0]->GetEntries();
2802 histopassedcuts[source*nbcuts+0]->SetStats(0);
2803 histopassedcuts[source*nbcuts+0]->Draw("same");
2804 leginput->AddEntry(histopassedcuts[source*nbcuts+0],"#pi^{0}","p");
2805 Double_t nbEta = 0.0;
2806 source = 3;
ccc37cdc 2807 if(((source*nbcuts+0)> (nbsources*nbcuts-1)) || (!histopassedcuts[source*nbcuts+0])) {
2808 delete [] histopassedcuts;
2809 delete [] nbEntriesCuts;
2810 return;
2811 }
e3ae862b 2812 nbEta = histopassedcuts[source*nbcuts+0]->GetEntries();
2813 histopassedcuts[source*nbcuts+0]->SetStats(0);
2814 histopassedcuts[source*nbcuts+0]->Draw("same");
2815 leginput->AddEntry(histopassedcuts[source*nbcuts+0],"#eta","p");
2816 Double_t nbC = 0.0;
2817 source = 4;
ccc37cdc 2818 if(((source*nbcuts+0)> (nbsources*nbcuts-1)) || (!histopassedcuts[source*nbcuts+0])) {
2819 delete [] histopassedcuts;
2820 delete [] nbEntriesCuts;
2821 return;
2822 }
e3ae862b 2823 nbC = histopassedcuts[source*nbcuts+0]->GetEntries();
2824 histopassedcuts[source*nbcuts+0]->SetStats(0);
2825 histopassedcuts[source*nbcuts+0]->Draw("same");
2826 leginput->AddEntry(histopassedcuts[source*nbcuts+0],"c","p");
2827 leginput->Draw("same");
2828
2829 //printf("Gamma %f, pi^{0} %f and #eta %f, c %f\n",nbGamma,nbPi0,nbEta,nbC);
2830
2831 //////////////////////
2832 // Tracked
2833 //////////////////////
2834
2835 TCanvas * cTracked = new TCanvas("cTracked","cTracked",800,800);
2836 cTracked->cd(1);
2837 source = 1;
ccc37cdc 2838 if(((source*nbcuts+1)> (nbsources*nbcuts-1)) || (!histopassedcuts[source*nbcuts+1])) {
2839 delete [] histopassedcuts;
2840 delete [] nbEntriesCuts;
2841 return;
2842 }
e3ae862b 2843 histopassedcuts[source*nbcuts+1]->Draw();
2844 TLegend *legTracked = new TLegend(0.4,0.6,0.89,0.89);
2845 legTracked->AddEntry(histopassedcuts[source*nbcuts+1],"#gamma","p");
2846 source = 2;
ccc37cdc 2847 if(((source*nbcuts+1)> (nbsources*nbcuts-1)) || (!histopassedcuts[source*nbcuts+1])) {
2848 delete [] histopassedcuts;
2849 delete [] nbEntriesCuts;
2850 return;
2851 }
e3ae862b 2852 histopassedcuts[source*nbcuts+1]->Draw("same");
2853 legTracked->AddEntry(histopassedcuts[source*nbcuts+1],"#pi^{0}","p");
2854 legTracked->Draw("same");
2855
ccc37cdc 2856 delete [] histopassedcuts;
2857 delete [] nbEntriesCuts;
e3ae862b 2858
70da6c5a 2859}
2860//_____________________________________________________________________________
2861Double_t AliHFEelecbackground::BetheBlochElectronITS(const Double_t *x, const Double_t * /*par*/)
2862{
2863 //
2864 // Bethe Bloch for ITS
2865 //
2866 static AliITSPIDResponse itsPidResponse;
2867 return itsPidResponse.Bethe(x[0],AliPID::ParticleMass(0));
2868}
2869//_____________________________________________________________________________
2870Double_t AliHFEelecbackground::BetheBlochMuonITS(const Double_t *x, const Double_t * /*par*/)
2871{
2872 //
2873 // Bethe Bloch for ITS
2874 //
2875 static AliITSPIDResponse itsPidResponse;
2876 return itsPidResponse.Bethe(x[0],AliPID::ParticleMass(1));
2877}
2878//_____________________________________________________________________________
2879Double_t AliHFEelecbackground::BetheBlochPionITS(const Double_t *x, const Double_t * /*par*/)
2880{
2881 //
2882 // Bethe Bloch for ITS
2883 //
2884 static AliITSPIDResponse itsPidResponse;
2885 return itsPidResponse.Bethe(x[0],AliPID::ParticleMass(2));
2886}
2887//_____________________________________________________________________________
2888Double_t AliHFEelecbackground::BetheBlochKaonITS(const Double_t *x, const Double_t * /*par*/)
2889{
2890 //
2891 // Bethe Bloch for ITS
2892 //
2893 static AliITSPIDResponse itsPidResponse;
2894 return itsPidResponse.Bethe(x[0],AliPID::ParticleMass(3));
2895}
2896//_____________________________________________________________________________
2897Double_t AliHFEelecbackground::BetheBlochProtonITS(const Double_t *x, const Double_t * /*par*/)
2898{
2899 //
2900 // Bethe Bloch for ITS
2901 //
2902 static AliITSPIDResponse itsPidResponse;
2903 return itsPidResponse.Bethe(x[0],AliPID::ParticleMass(4));
2904}
2905//_____________________________________________________________________________
2906Double_t AliHFEelecbackground::BetheBlochElectronTPC(const Double_t *x, const Double_t * /*par*/)
2907{
2908 //
2909 // Bethe Bloch for TPC
2910 //
2911 static AliTPCPIDResponse tpcPidResponse;
2912 return tpcPidResponse.GetExpectedSignal(x[0],AliPID::kElectron);
2913}
2914//_____________________________________________________________________________
2915Double_t AliHFEelecbackground::BetheBlochMuonTPC(const Double_t *x, const Double_t * /*par*/)
2916{
2917 //
2918 // Bethe Bloch for TPC
2919 //
2920 static AliTPCPIDResponse tpcPidResponse;
2921 return tpcPidResponse.GetExpectedSignal(x[0],AliPID::kMuon);
2922}
2923//_____________________________________________________________________________
2924Double_t AliHFEelecbackground::BetheBlochPionTPC(const Double_t *x, const Double_t * /*par*/)
2925{
2926 //
2927 // Bethe Bloch for TPC
2928 //
2929 static AliTPCPIDResponse tpcPidResponse;
2930 return tpcPidResponse.GetExpectedSignal(x[0],AliPID::kPion);
2931}
2932//_____________________________________________________________________________
2933Double_t AliHFEelecbackground::BetheBlochKaonTPC(const Double_t *x, const Double_t * /*par*/)
2934{
2935 //
2936 // Bethe Bloch for TPC
2937 //
2938 static AliTPCPIDResponse tpcPidResponse;
2939 return tpcPidResponse.GetExpectedSignal(x[0],AliPID::kKaon);
2940}
2941//_____________________________________________________________________________
2942Double_t AliHFEelecbackground::BetheBlochProtonTPC(const Double_t *x, const Double_t * /*par*/)
2943{
2944 //
2945 // Bethe Bloch for TPC
2946 //
2947 static AliTPCPIDResponse tpcPidResponse;
2948 return tpcPidResponse.GetExpectedSignal(x[0],AliPID::kProton);
02524e30 2949}
70da6c5a 2950
2951
2952