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