New class for normalization studies (Giacomo)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAnalysisTaskSED0Mass.cxx
CommitLineData
49061176 1/**************************************************************************
2 * Copyright(c) 1998-2009, 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//
18// AliAnalysisTaskSE for D0 candidates invariant mass histogram
a41f6fad 19// and comparison with the MC truth and cut variables distributions.
49061176 20//
21// Authors: A.Dainese, andrea.dainese@lnl.infn.it
feb73eca 22// Chiara Bianchin, chiara.bianchin@pd.infn.it (invariant mass)
23// Carmelo Di Giglio, carmelo.digiglio@ba.infn.it (like sign)
49061176 24/////////////////////////////////////////////////////////////
25
26#include <Riostream.h>
27#include <TClonesArray.h>
527f330b 28#include <TCanvas.h>
49061176 29#include <TNtuple.h>
30#include <TList.h>
31#include <TH1F.h>
a41f6fad 32#include <TH2F.h>
33#include <TDatabasePDG.h>
49061176 34
5b2e5fae 35#include <AliAnalysisDataSlot.h>
36#include <AliAnalysisDataContainer.h>
b557eb43 37#include "AliAnalysisManager.h"
34137226 38#include "AliESDtrack.h"
4e61a020 39#include "AliVertexerTracks.h"
b557eb43 40#include "AliAODHandler.h"
49061176 41#include "AliAODEvent.h"
42#include "AliAODVertex.h"
43#include "AliAODTrack.h"
44#include "AliAODMCHeader.h"
45#include "AliAODMCParticle.h"
46#include "AliAODRecoDecayHF2Prong.h"
47#include "AliAODRecoCascadeHF.h"
48#include "AliAnalysisVertexingHF.h"
49#include "AliAnalysisTaskSE.h"
50#include "AliAnalysisTaskSED0Mass.h"
a96083b9 51#include "AliNormalizationCounter.h"
b557eb43 52
49061176 53ClassImp(AliAnalysisTaskSED0Mass)
54
55
56//________________________________________________________________________
57AliAnalysisTaskSED0Mass::AliAnalysisTaskSED0Mass():
58AliAnalysisTaskSE(),
ea0d8716 59fOutputMass(0),
60fDistr(0),
61fNentries(0),
34137226 62fChecks(0),
ea0d8716 63fCuts(0),
feb73eca 64fArray(0),
ce39f0ac 65fReadMC(0),
40445ada 66fCutOnDistr(0),
6b3e3c78 67fUsePid4Distr(0),
a96083b9 68fCounter(0),
87020237 69fNPtBins(1),
4b4a1d25 70fTotPosPairs(0),
71fTotNegPairs(0),
449b1302 72fLsNormalization(1.),
73fFillOnlyD0D0bar(0)
feb73eca 74
ea0d8716 75
49061176 76{
77 // Default constructor
78}
79
80//________________________________________________________________________
ea0d8716 81AliAnalysisTaskSED0Mass::AliAnalysisTaskSED0Mass(const char *name,AliRDHFCutsD0toKpi* cuts):
49061176 82AliAnalysisTaskSE(name),
ea0d8716 83fOutputMass(0),
34137226 84fDistr(0),
ea0d8716 85fNentries(0),
34137226 86fChecks(0),
ea0d8716 87fCuts(0),
feb73eca 88fArray(0),
ce39f0ac 89fReadMC(0),
40445ada 90fCutOnDistr(0),
6b3e3c78 91fUsePid4Distr(0),
a96083b9 92fCounter(0),
87020237 93fNPtBins(1),
4b4a1d25 94fTotPosPairs(0),
95fTotNegPairs(0),
449b1302 96fLsNormalization(1.),
97fFillOnlyD0D0bar(0)
ea0d8716 98
49061176 99{
100 // Default constructor
87020237 101
102 fNPtBins=cuts->GetNPtBins();
103 fTotPosPairs=new Int_t[fNPtBins];
104 fTotNegPairs=new Int_t[fNPtBins];
105 for(Int_t i=0;i<fNPtBins;i++) {fTotPosPairs[i]=0; fTotNegPairs[i]=0;}
4b4a1d25 106
ea0d8716 107 fCuts=cuts;
108
109 // Output slot #1 writes into a TList container (mass with cuts)
49061176 110 DefineOutput(1,TList::Class()); //My private output
ea0d8716 111 // Output slot #2 writes into a TList container (distributions)
a4ae02cd 112 DefineOutput(2,TList::Class()); //My private output
ea0d8716 113 // Output slot #3 writes into a TH1F container (number of events)
a4ae02cd 114 DefineOutput(3,TH1F::Class()); //My private output
ea0d8716 115 // Output slot #4 writes into a TList container (quality check)
a41f6fad 116 DefineOutput(4,TList::Class()); //My private output
ea0d8716 117 // Output slot #5 writes into a TList container (cuts)
118 DefineOutput(5,AliRDHFCutsD0toKpi::Class()); //My private output
a96083b9 119 // Output slot #6 writes Normalization Counter
120 DefineOutput(6,AliNormalizationCounter::Class());
49061176 121}
122
123//________________________________________________________________________
124AliAnalysisTaskSED0Mass::~AliAnalysisTaskSED0Mass()
125{
ea0d8716 126 if (fOutputMass) {
127 delete fOutputMass;
128 fOutputMass = 0;
a4ae02cd 129 }
a41f6fad 130 if (fDistr) {
131 delete fDistr;
132 fDistr = 0;
133 }
34137226 134 if (fChecks) {
135 delete fChecks;
136 fChecks = 0;
137 }
ea0d8716 138 if (fCuts) {
139 delete fCuts;
140 fCuts = 0;
49061176 141 }
46c96ce5 142 if (fNentries){
143 delete fNentries;
144 fNentries = 0;
145 }
a96083b9 146 if(fCounter){
147 delete fCounter;
148 fCounter=0;
149 }
ea0d8716 150
49061176 151}
152
153//________________________________________________________________________
154void AliAnalysisTaskSED0Mass::Init()
155{
156 // Initialization
157
158 if(fDebug > 1) printf("AnalysisTaskSED0Mass::Init() \n");
159
ea0d8716 160
ea0d8716 161 AliRDHFCutsD0toKpi* copyfCuts=new AliRDHFCutsD0toKpi(*fCuts);
4e61a020 162 const char* nameoutput=GetOutputSlot(5)->GetContainer()->GetName();
163 copyfCuts->SetName(nameoutput);
ea0d8716 164 // Post the data
165 PostData(5,copyfCuts);
4e61a020 166
49061176 167
49061176 168 return;
169}
170
171//________________________________________________________________________
172void AliAnalysisTaskSED0Mass::UserCreateOutputObjects()
173{
174
175 // Create the output container
176 //
177 if(fDebug > 1) printf("AnalysisTaskSED0Mass::UserCreateOutputObjects() \n");
178
179 // Several histograms are more conveniently managed in a TList
ea0d8716 180 fOutputMass = new TList();
181 fOutputMass->SetOwner();
182 fOutputMass->SetName("listMass");
49061176 183
a41f6fad 184 fDistr = new TList();
185 fDistr->SetOwner();
186 fDistr->SetName("distributionslist");
187
34137226 188 fChecks = new TList();
189 fChecks->SetOwner();
190 fChecks->SetName("checkHistograms");
191
0108fa62 192 TString nameMass=" ",nameSgn27=" ",nameSgn=" ", nameBkg=" ", nameRfl=" ",nameMassNocutsS =" ",nameMassNocutsB =" ", namedistr=" ";
49061176 193
ea0d8716 194 for(Int_t i=0;i<fCuts->GetNPtBins();i++){
b272aebf 195
7646d6da 196 nameMass="histMass_";
ea0d8716 197 nameMass+=i;
9de8c723 198 nameSgn27="histSgn27_";
ea0d8716 199 nameSgn27+=i;
7646d6da 200 nameSgn="histSgn_";
ea0d8716 201 nameSgn+=i;
7646d6da 202 nameBkg="histBkg_";
ea0d8716 203 nameBkg+=i;
a4ae02cd 204 nameRfl="histRfl_";
ea0d8716 205 nameRfl+=i;
0108fa62 206 nameMassNocutsS="hMassS_";
ea0d8716 207 nameMassNocutsS+=i;
0108fa62 208 nameMassNocutsB="hMassB_";
ea0d8716 209 nameMassNocutsB+=i;
7646d6da 210
b272aebf 211 //histograms of cut variable distributions
212
213 // pT
214 namedistr="hptpiS_";
ea0d8716 215 namedistr+=i;
b272aebf 216 TH1F *hptpiS = new TH1F(namedistr.Data(), "P_{T} distribution (pions);p_{T} [GeV/c]",200,0.,8.);
217
218 namedistr="hptKS_";
ea0d8716 219 namedistr+=i;
b272aebf 220 TH1F *hptKS = new TH1F(namedistr.Data(), "P_{T} distribution (kaons);p_{T} [GeV/c]",200,0.,8.);
221
222 namedistr="hptB_";
ea0d8716 223 namedistr+=i;
b272aebf 224 TH1F *hptB = new TH1F(namedistr.Data(), "P_{T} distribution;p_{T} [GeV/c]",200,0.,8.);
225
226 // pT no mass cut
227 namedistr="hptpiSnoMcut_";
ea0d8716 228 namedistr+=i;
b272aebf 229 TH1F *hptpiSnoMcut = new TH1F(namedistr.Data(), "P_{T} distribution (pions);p_{T} [GeV/c]",200,0.,8.);
230
231 namedistr="hptKSnoMcut_";
ea0d8716 232 namedistr+=i;
b272aebf 233 TH1F *hptKSnoMcut = new TH1F(namedistr.Data(), "P_{T} distribution (kaons);p_{T} [GeV/c]",200,0.,8.);
234
235 namedistr="hptB1prongnoMcut_";
ea0d8716 236 namedistr+=i;
b272aebf 237 TH1F *hptB1pnoMcut = new TH1F(namedistr.Data(), "P_{T} distribution;p_{T} [GeV/c]",200,0.,8.);
238
239 namedistr="hptB2prongsnoMcut_";
ea0d8716 240 namedistr+=i;
b272aebf 241 TH1F *hptB2pnoMcut = new TH1F(namedistr.Data(), "P_{T} distribution;p_{T} [GeV/c]",200,0.,8.);
242
243
244 // dca
245 namedistr="hdcaS_";
ea0d8716 246 namedistr+=i;
b272aebf 247 TH1F *hdcaS = new TH1F(namedistr.Data(), "DCA distribution;dca [cm]",200,0.,0.1);
248 namedistr="hdcaB_";
ea0d8716 249 namedistr+=i;
b272aebf 250 TH1F *hdcaB = new TH1F(namedistr.Data(), "DCA distribution;dca [cm]",200,0.,0.1);
251
252 // costhetastar
253 namedistr="hcosthetastarS_";
ea0d8716 254 namedistr+=i;
b272aebf 255 TH1F *hcosthetastarS = new TH1F(namedistr.Data(), "cos#theta* distribution;cos#theta*",200,-1.,1.);
256 namedistr="hcosthetastarB_";
ea0d8716 257 namedistr+=i;
b272aebf 258 TH1F *hcosthetastarB = new TH1F(namedistr.Data(), "cos#theta* distribution;cos#theta*",200,-1.,1.);
259
260 // impact parameter
261 namedistr="hd0piS_";
ea0d8716 262 namedistr+=i;
b272aebf 263 TH1F *hd0piS = new TH1F(namedistr.Data(), "Impact parameter distribution (pions);d0(#pi) [cm]",200,-0.1,0.1);
264
265 namedistr="hd0KS_";
ea0d8716 266 namedistr+=i;
b272aebf 267 TH1F *hd0KS = new TH1F(namedistr.Data(), "Impact parameter distribution (kaons);d0(K) [cm]",200,-0.1,0.1);
4e61a020 268
b272aebf 269 namedistr="hd0B_";
ea0d8716 270 namedistr+=i;
4e61a020 271 TH1F *hd0B = new TH1F(namedistr.Data(), "Impact parameter distribution (both);d0 [cm]",200,-0.1,0.1);
272
273 namedistr="hd0p0B_";
274 namedistr+=i;
275 TH1F *hd0p0B = new TH1F(namedistr.Data(), "Impact parameter distribution (prong +);d0 [cm]",200,-0.1,0.1);
276
6b3e3c78 277 namedistr="hd0moresB_";
278 namedistr+=i;
279 TH1F *hd0moresB = new TH1F(namedistr.Data(), "Impact parameter distribution (both);d0 [cm]",200,-0.1,0.1);
280
281 namedistr="hd0d0moresB_";
282 namedistr+=i;
283 TH1F *hd0d0moresB = new TH1F(namedistr.Data(), "Impact parameter distribution (prong +);d0 [cm]",200,-0.001,0.001);
284
285 namedistr="hd0vmoresB_";
286 namedistr+=i;
287 TH1F *hd0vmoresB = new TH1F(namedistr.Data(), "Impact parameter distribution (both);d0 [cm]",200,-0.1,0.1);
288
289 namedistr="hd0d0vmoresB_";
290 namedistr+=i;
291 TH1F *hd0d0vmoresB = new TH1F(namedistr.Data(), "Impact parameter distribution (prong +);d0 [cm]",200,-0.001,0.001);
292
4e61a020 293 namedistr="hd0p1B_";
294 namedistr+=i;
295 TH1F *hd0p1B = new TH1F(namedistr.Data(), "Impact parameter distribution (prong -);d0 [cm]",200,-0.1,0.1);
296
297 namedistr="hd0vpiS_";
298 namedistr+=i;
299 TH1F *hd0vpiS = new TH1F(namedistr.Data(), "Impact parameter distribution (pions)(vtx w/o these tracks);d0(#pi) [cm]",200,-0.1,0.1);
300
301 namedistr="hd0vKS_";
302 namedistr+=i;
303 TH1F *hd0vKS = new TH1F(namedistr.Data(), "Impact parameter distribution (kaons) (vtx w/o these tracks);d0(K) [cm]",200,-0.1,0.1);
304 namedistr="hd0vB_";
305 namedistr+=i;
306 TH1F *hd0vB = new TH1F(namedistr.Data(), "Impact parameter distribution (vtx w/o these tracks);d0 [cm]",200,-0.1,0.1);
307
308 namedistr="hd0vp0B_";
309 namedistr+=i;
310 TH1F *hd0vp0B = new TH1F(namedistr.Data(), "Impact parameter distribution (prong + ** vtx w/o these tracks);d0 [cm]",200,-0.1,0.1);
311
312 namedistr="hd0vp1B_";
313 namedistr+=i;
314 TH1F *hd0vp1B = new TH1F(namedistr.Data(), "Impact parameter distribution (prong - ** vtx w/o these tracks);d0 [cm]",200,-0.1,0.1);
b272aebf 315
316 namedistr="hd0d0S_";
ea0d8716 317 namedistr+=i;
b272aebf 318 TH1F *hd0d0S = new TH1F(namedistr.Data(), "d_{0}#timesd_{0} distribution;d_{0}#timesd_{0} [cm^{2}]",200,-0.001,0.001);
319 namedistr="hd0d0B_";
ea0d8716 320 namedistr+=i;
b272aebf 321 TH1F *hd0d0B = new TH1F(namedistr.Data(), "d_{0}#timesd_{0} distribution;d_{0}#timesd_{0} [cm^{2}]",200,-0.001,0.001);
322
4e61a020 323 namedistr="hd0d0vS_";
324 namedistr+=i;
325 TH1F *hd0d0vS = new TH1F(namedistr.Data(), "d_{0}#timesd_{0} distribution (vtx w/o these tracks);d_{0}#timesd_{0} [cm^{2}]",200,-0.001,0.001);
326 namedistr="hd0d0vB_";
327 namedistr+=i;
328 TH1F *hd0d0vB = new TH1F(namedistr.Data(), "d_{0}#timesd_{0} distribution (vtx w/o these tracks);d_{0}#timesd_{0} [cm^{2}]",200,-0.001,0.001);
329
330 //decay lenght
331 namedistr="hdeclS_";
332 namedistr+=i;
6b3e3c78 333 TH1F *hdeclengthS = new TH1F(namedistr.Data(), "Decay Length distribution;Decay Length [cm]",200,0,0.15);
4e61a020 334
335 namedistr="hdeclB_";
336 namedistr+=i;
6b3e3c78 337 TH1F *hdeclengthB = new TH1F(namedistr.Data(), "Decay Length distribution;Decay Length [cm]",200,0,0.15);
4e61a020 338
339 namedistr="hnormdeclS_";
340 namedistr+=i;
6b3e3c78 341 TH1F *hnormdeclengthS = new TH1F(namedistr.Data(), "Normalized Decay Length distribution;Decay Length/Err ",200,0,6.);
4e61a020 342
343 namedistr="hnormdeclB_";
344 namedistr+=i;
6b3e3c78 345 TH1F *hnormdeclengthB = new TH1F(namedistr.Data(), "Normalized Decay Length distribution;Decay Length/Err ",200,0,6.);
4e61a020 346
a2121012 347 namedistr="hdeclvS_";
348 namedistr+=i;
349 TH1F *hdeclengthvS = new TH1F(namedistr.Data(), "Decay Length distribution (vtx w/o tracks);Decay Length [cm]",200,0,0.6);
350
351 namedistr="hdeclvB_";
352 namedistr+=i;
353 TH1F *hdeclengthvB = new TH1F(namedistr.Data(), "Decay Length distribution (vtx w/o tracks);Decay Length [cm]",200,0,0.6);
354
355 namedistr="hnormdeclvS_";
356 namedistr+=i;
357 TH1F *hnormdeclengthvS = new TH1F(namedistr.Data(), "Normalized Decay Length distribution (vtx w/o tracks);Decay Length/Err ",200,0,10.);
358
359 namedistr="hnormdeclvB_";
360 namedistr+=i;
361 TH1F *hnormdeclengthvB = new TH1F(namedistr.Data(), "Normalized Decay Length distribution (vtx w/o tracks);Decay Length/Err ",200,0,10.);
362
4e61a020 363 // costhetapoint
b272aebf 364 namedistr="hcosthetapointS_";
ea0d8716 365 namedistr+=i;
b272aebf 366 TH1F *hcosthetapointS = new TH1F(namedistr.Data(), "cos#theta_{Point} distribution;cos#theta_{Point}",200,0,1.);
367 namedistr="hcosthetapointB_";
ea0d8716 368 namedistr+=i;
b272aebf 369 TH1F *hcosthetapointB = new TH1F(namedistr.Data(), "cos#theta_{Point} distribution;cos#theta_{Point}",200,0,1.);
370
6b3e3c78 371 namedistr="hcosthetapointmoresB_";
372 namedistr+=i;
373 TH1F *hcosthetapointmoresB = new TH1F(namedistr.Data(), "cos#theta_{Point} distribution;cos#theta_{Point}",200,0,1.);
374
4e61a020 375 // costhetapoint vs d0 or d0d0
376 namedistr="hcosthpointd0S_";
377 namedistr+=i;
378 TH2F *hcosthpointd0S= new TH2F(namedistr.Data(),"Correlation cos#theta_{Point}-d_{0};cos#theta_{Point};d_{0} [cm^{2}]",200,0,1.,200,-0.001,0.001);
379
380 namedistr="hcosthpointd0B_";
381 namedistr+=i;
382 TH2F *hcosthpointd0B= new TH2F(namedistr.Data(),"Correlation cos#theta_{Point}-d_{0};cos#theta_{Point};d_{0} [cm^{2}]",200,0,1.,200,-0.001,0.001);
383
b272aebf 384 namedistr="hcosthpointd0d0S_";
ea0d8716 385 namedistr+=i;
b272aebf 386 TH2F *hcosthpointd0d0S= new TH2F(namedistr.Data(),"Correlation cos#theta_{Point}-d_{0}#timesd_{0};cos#theta_{Point};d_{0}#timesd_{0} [cm^{2}]",200,0,1.,200,-0.001,0.001);
387 namedistr="hcosthpointd0d0B_";
ea0d8716 388 namedistr+=i;
b272aebf 389 TH2F *hcosthpointd0d0B= new TH2F(namedistr.Data(),"Correlation cos#theta_{Point}-d_{0}#timesd_{0};cos#theta_{Point};d_{0}#timesd_{0} [cm^{2}]",200,0,1.,200,-0.001,0.001);
390
391 fDistr->Add(hptpiS);
392 fDistr->Add(hptKS);
393 fDistr->Add(hptB);
394
395 fDistr->Add(hptpiSnoMcut);
396 fDistr->Add(hptKSnoMcut);
397 fDistr->Add(hptB1pnoMcut);
398 fDistr->Add(hptB2pnoMcut);
399
400 fDistr->Add(hdcaS);
401 fDistr->Add(hdcaB);
402
403 fDistr->Add(hd0piS);
404 fDistr->Add(hd0KS);
405 fDistr->Add(hd0B);
4e61a020 406 fDistr->Add(hd0p0B);
6b3e3c78 407 fDistr->Add(hd0moresB);
4e61a020 408 fDistr->Add(hd0p1B);
409
410 fDistr->Add(hd0vpiS);
411 fDistr->Add(hd0vKS);
412 fDistr->Add(hd0vB);
413 fDistr->Add(hd0vp0B);
414 fDistr->Add(hd0vp1B);
6b3e3c78 415 fDistr->Add(hd0vmoresB);
b272aebf 416
417 fDistr->Add(hd0d0S);
418 fDistr->Add(hd0d0B);
6b3e3c78 419 fDistr->Add(hd0d0moresB);
b272aebf 420
4e61a020 421 fDistr->Add(hd0d0vS);
422 fDistr->Add(hd0d0vB);
6b3e3c78 423 fDistr->Add(hd0d0vmoresB);
4e61a020 424
b272aebf 425 fDistr->Add(hcosthetastarS);
426 fDistr->Add(hcosthetastarB);
427
428 fDistr->Add(hcosthetapointS);
429 fDistr->Add(hcosthetapointB);
6b3e3c78 430 fDistr->Add(hcosthetapointmoresB);
b272aebf 431
4e61a020 432 fDistr->Add(hdeclengthS);
433 fDistr->Add(hdeclengthB);
434
435 fDistr->Add(hnormdeclengthS);
436 fDistr->Add(hnormdeclengthB);
437
a2121012 438 fDistr->Add(hdeclengthvS);
439 fDistr->Add(hdeclengthvB);
440
441 fDistr->Add(hnormdeclengthvS);
442 fDistr->Add(hnormdeclengthvB);
443
4e61a020 444 fDistr->Add(hcosthpointd0S);
445 fDistr->Add(hcosthpointd0B);
446
b272aebf 447 fDistr->Add(hcosthpointd0d0S);
448 fDistr->Add(hcosthpointd0d0B);
449
450
a41f6fad 451 //histograms of invariant mass distributions
452
46c96ce5 453 TH1F* tmpMt = new TH1F(nameMass.Data(),"D^{0} invariant mass; M [GeV]; Entries",200,1.765,1.965);
454 TH1F *tmpMl=(TH1F*)tmpMt->Clone();
455 tmpMt->Sumw2();
456 tmpMl->Sumw2();
a4ae02cd 457
9de8c723 458 //to compare with AliAnalysisTaskCharmFraction
459 TH1F* tmpS27t = new TH1F(nameSgn27.Data(),"D^{0} invariant mass in M(D^{0}) +/- 27 MeV - MC; M [GeV]; Entries",200,1.765,1.965);
460 TH1F *tmpS27l=(TH1F*)tmpS27t->Clone();
461 tmpS27t->Sumw2();
462 tmpS27l->Sumw2();
0108fa62 463
464 //distribution w/o cuts
6321ee46 465 // TH1F* tmpMS = new TH1F(nameMassNocutsS.Data(),"D^{0} invariant mass; M [GeV]; Entries",300,0.7,3.);
466 TH1F* tmpMS = new TH1F(nameMassNocutsS.Data(),"D^{0} invariant mass; M [GeV]; Entries",300,1.56484,2.16484); //range (MD0-300MeV, mD0 + 300MeV)
467 TH1F *tmpMB=(TH1F*)tmpMS->Clone();
0108fa62 468 tmpMB->SetName(nameMassNocutsB.Data());
469 tmpMS->Sumw2();
470 tmpMB->Sumw2();
471
472 //MC signal and background
46c96ce5 473 TH1F* tmpSt = new TH1F(nameSgn.Data(), "D^{0} invariant mass - MC; M [GeV]; Entries",200,1.765,1.965);
474 TH1F *tmpSl=(TH1F*)tmpSt->Clone();
475 tmpSt->Sumw2();
476 tmpSl->Sumw2();
a4ae02cd 477
46c96ce5 478 TH1F* tmpBt = new TH1F(nameBkg.Data(), "Background invariant mass - MC; M [GeV]; Entries",200,1.765,1.965);
479 TH1F *tmpBl=(TH1F*)tmpBt->Clone();
480 tmpBt->Sumw2();
481 tmpBl->Sumw2();
a4ae02cd 482
feb73eca 483 //Reflection: histo filled with D0Mass which pass the cut (also) as D0bar and with D0bar which pass (also) the cut as D0
46c96ce5 484 TH1F* tmpRt = new TH1F(nameRfl.Data(), "Reflected signal invariant mass - MC; M [GeV]; Entries",200,1.765,1.965);
485 TH1F *tmpRl=(TH1F*)tmpRt->Clone();
486 tmpRt->Sumw2();
487 tmpRl->Sumw2();
a4ae02cd 488 // printf("Created histograms %s\t%s\t%s\t%s\n",tmpM->GetName(),tmpS->GetName(),tmpB->GetName(),tmpR->GetName());
489
ea0d8716 490 fOutputMass->Add(tmpMt);
491 fOutputMass->Add(tmpSt);
492 fOutputMass->Add(tmpS27t);
493 fOutputMass->Add(tmpBt);
494 fOutputMass->Add(tmpRt);
495
0108fa62 496 fDistr->Add(tmpMS);
497 fDistr->Add(tmpMB);
498
34137226 499
7646d6da 500 }
a4ae02cd 501
449b1302 502 namedistr="hpospair";
503 TH1F* hpospair=new TH1F(namedistr.Data(),"Number of positive pairs",fCuts->GetNPtBins(),-0.5,fCuts->GetNPtBins()-0.5);
504 namedistr="hnegpair";
505 TH1F* hnegpair=new TH1F(namedistr.Data(),"Number of negative pairs",fCuts->GetNPtBins(),-0.5,fCuts->GetNPtBins()-0.5);
506 fDistr->Add(hpospair);
507 fDistr->Add(hnegpair);
508
509
4e61a020 510 //histograms for vertex checking and TOF checking
34137226 511 TString checkname="hptGoodTr";
512
513 TH1F* hptGoodTr=new TH1F(checkname.Data(),"Pt distribution of 'good' tracks;p_{t}[GeV];Number",200,0.,8.);
514 hptGoodTr->SetTitleOffset(1.3,"Y");
515 checkname="hdistrGoodTr";
a41f6fad 516
34137226 517 TH1F* hdistrGoodTr=new TH1F(checkname.Data(),"Distribution of number of good tracks per event;no.good-tracks/ev;Entries",31,0,31);
518 hdistrGoodTr->SetTitleOffset(1.3,"Y");
519
4e61a020 520 checkname="hTOFsig";
521 TH1F* hTOFsig=new TH1F(checkname.Data(),"Distribution of TOF signal;TOF time [ps];Entries", 100, -2.e3,40.e3);
522
6b3e3c78 523 checkname="hTPCsig";
524 TH1F* hTPCsig=new TH1F(checkname.Data(),"Distribution of TPC signal;TPC sig;Entries", 100, 35.,100.);
525
4e61a020 526 checkname="hTOFtimeKaonHyptime";
6b3e3c78 527 TH2F* hTOFtimeKaonHyptime=new TH2F(checkname.Data(),"TOFtime - timeHypothesisForKaon;p[GeV/c];TOFtime - timeHypothesisForKaon [ps]",200,0.,4.,1000,-20000.,20000.);
528
529 checkname="hTOFtimeKaonHyptimeAC";
530 TH2F* hTOFtimeKaonHyptimeAC=new TH2F(checkname.Data(),"TOFtime - timeHypothesisForKaon;p[GeV/c];TOFtime - timeHypothesisForKaon [ps]",200,0.,4.,1000,-20000.,20000.);
531
532 checkname="hTPCsigvsp";
533 TH2F* hTPCsigvsp=new TH2F(checkname.Data(),"TPCsig vs p;p[GeV/c];TPCsig",200,0.,4.,1000,35.,100.);
4e61a020 534
535 checkname="hTOFtime";
536 TH1F* hTOFtime=new TH1F(checkname.Data(),"Distribution of TOF time Kaon;TOF time(Kaon) [ps];Entries", 1000, 0.,50000.);
537
538
34137226 539 fChecks->Add(hptGoodTr);
540 fChecks->Add(hdistrGoodTr);
4e61a020 541 fChecks->Add(hTOFsig);
6b3e3c78 542 fChecks->Add(hTPCsig);
4e61a020 543 fChecks->Add(hTOFtimeKaonHyptime);
6b3e3c78 544 fChecks->Add(hTOFtimeKaonHyptimeAC);
545 fChecks->Add(hTPCsigvsp);
4e61a020 546 fChecks->Add(hTOFtime);
34137226 547
5b2e5fae 548 const char* nameoutput=GetOutputSlot(3)->GetContainer()->GetName();
4e61a020 549
6b3e3c78 550 fNentries=new TH1F(nameoutput, "Integral(1,2) = number of AODs *** Integral(2,3) = number of candidates selected with cuts *** Integral(3,4) = number of D0 selected with cuts *** Integral(4,5) = events with good vertex *** Integral(5,6) = pt out of bounds", 14,-0.5,13.5);
34137226 551
552 fNentries->GetXaxis()->SetBinLabel(1,"nEventsAnal");
4e61a020 553 fNentries->GetXaxis()->SetBinLabel(2,"nCandSel(Cuts)");
34137226 554 fNentries->GetXaxis()->SetBinLabel(3,"nD0Selected");
6b3e3c78 555 fNentries->GetXaxis()->SetBinLabel(4,"nEventsGoodVtxS");
34137226 556 fNentries->GetXaxis()->SetBinLabel(5,"nEventsGoodVtx+>2tracks");
4e61a020 557 fNentries->GetXaxis()->SetBinLabel(6,"ptbin = -1");
558 fNentries->GetXaxis()->SetBinLabel(7,"no daughter");
559 fNentries->GetXaxis()->SetBinLabel(8,"nCandSel(Tr)");
6b3e3c78 560 fNentries->GetXaxis()->SetBinLabel(9,"PID=0");
561 fNentries->GetXaxis()->SetBinLabel(10,"PID=1");
562 fNentries->GetXaxis()->SetBinLabel(11,"PID=2");
563 fNentries->GetXaxis()->SetBinLabel(12,"PID=3");
564 fNentries->GetXaxis()->SetBinLabel(13,"K");
565 fNentries->GetXaxis()->SetBinLabel(14,"Lambda");
34137226 566 fNentries->GetXaxis()->SetNdivisions(1,kFALSE);
a4ae02cd 567
a96083b9 568 fCounter = new AliNormalizationCounter("NormalizationCounter");
569
40445ada 570 // Post the data
ea0d8716 571 PostData(1,fOutputMass);
5b2e5fae 572 PostData(2,fDistr);
40445ada 573 PostData(3,fNentries);
5b2e5fae 574 PostData(4,fChecks);
a96083b9 575 PostData(6,fCounter);
49061176 576 return;
577}
578
579//________________________________________________________________________
580void AliAnalysisTaskSED0Mass::UserExec(Option_t */*option*/)
581{
582 // Execute analysis for current event:
583 // heavy flavor candidates association to MC truth
a4ae02cd 584 //cout<<"I'm in UserExec"<<endl;
ea0d8716 585
586
587 //cuts order
588 // printf(" |M-MD0| [GeV] < %f\n",fD0toKpiCuts[0]);
589 // printf(" dca [cm] < %f\n",fD0toKpiCuts[1]);
590 // printf(" cosThetaStar < %f\n",fD0toKpiCuts[2]);
591 // printf(" pTK [GeV/c] > %f\n",fD0toKpiCuts[3]);
592 // printf(" pTpi [GeV/c] > %f\n",fD0toKpiCuts[4]);
593 // printf(" |d0K| [cm] < %f\n",fD0toKpiCuts[5]);
594 // printf(" |d0pi| [cm] < %f\n",fD0toKpiCuts[6]);
595 // printf(" d0d0 [cm^2] < %f\n",fD0toKpiCuts[7]);
596 // printf(" cosThetaPoint > %f\n",fD0toKpiCuts[8]);
597
598
49061176 599 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
feb73eca 600
b557eb43 601 TString bname;
feb73eca 602 if(fArray==0){ //D0 candidates
b557eb43 603 // load D0->Kpi candidates
feb73eca 604 //cout<<"D0 candidates"<<endl;
b557eb43 605 bname="D0toKpi";
feb73eca 606 } else { //LikeSign candidates
feb73eca 607 // load like sign candidates
b557eb43 608 //cout<<"LS candidates"<<endl;
609 bname="LikeSign2Prong";
610 }
611
612 TClonesArray *inputArray=0;
34137226 613
b557eb43 614 if(!aod && AODEvent() && IsStandardAOD()) {
615 // In case there is an AOD handler writing a standard AOD, use the AOD
616 // event in memory rather than the input (ESD) event.
617 aod = dynamic_cast<AliAODEvent*> (AODEvent());
618 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
619 // have to taken from the AOD event hold by the AliAODExtension
620 AliAODHandler* aodHandler = (AliAODHandler*)
621 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
34137226 622
b557eb43 623 if(aodHandler->GetExtensions()) {
624 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
34137226 625 AliAODEvent* aodFromExt = ext->GetAOD();
b557eb43 626 inputArray=(TClonesArray*)aodFromExt->GetList()->FindObject(bname.Data());
feb73eca 627 }
b557eb43 628 } else {
629 inputArray=(TClonesArray*)aod->GetList()->FindObject(bname.Data());
630 }
feb73eca 631
b557eb43 632
633 if(!inputArray) {
634 printf("AliAnalysisTaskSED0Mass::UserExec: input branch not found!\n");
635 return;
49061176 636 }
637
6b3e3c78 638 // fix for temporary bug in ESDfilter
7c23877d 639 // the AODs with null vertex pointer didn't pass the PhysSel
5806c290 640 if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) return;
641
7c23877d 642
ce39f0ac 643 TClonesArray *mcArray = 0;
644 AliAODMCHeader *mcHeader = 0;
645
646 if(fReadMC) {
647 // load MC particles
648 mcArray = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
649 if(!mcArray) {
650 printf("AliAnalysisTaskSED0Mass::UserExec: MC particles branch not found!\n");
651 return;
652 }
40445ada 653
ce39f0ac 654 // load MC header
655 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
656 if(!mcHeader) {
657 printf("AliAnalysisTaskSED0Mass::UserExec: MC header branch not found!\n");
658 return;
659 }
49061176 660 }
661
b557eb43 662 //printf("VERTEX Z %f %f\n",vtx1->GetZ(),mcHeader->GetVtxZ());
40445ada 663
a4ae02cd 664 //histogram filled with 1 for every AOD
34137226 665 fNentries->Fill(0);
a96083b9 666 fCounter->StoreEvent(aod,fReadMC);
6b3e3c78 667 if(!fCuts->IsEventSelected(aod)) return;
40445ada 668
669 // AOD primary vertex
670 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
b272aebf 671
40445ada 672 Bool_t isGoodVtx=kFALSE;
b272aebf 673
40445ada 674 //vtx1->Print();
675 TString primTitle = vtx1->GetTitle();
676 if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0) {
677 isGoodVtx=kTRUE;
678 fNentries->Fill(3);
679 }
680
6b3e3c78 681
40445ada 682 //cout<<"Start checks"<<endl;
683 Int_t ntracks=0,isGoodTrack=0;
684
685 if(aod) ntracks=aod->GetNTracks();
686
75638da0 687 //cout<<"ntracks = "<<ntracks<<endl;
ea0d8716 688
75638da0 689 //loop on tracks in the event
690 for (Int_t k=0;k<ntracks;k++){
691 AliAODTrack* track=aod->GetTrack(k);
6b3e3c78 692 //cout<<"track n = "<<k<<endl;
693
4e61a020 694 //check TOF
4e61a020 695 if(!(track->GetStatus()&AliESDtrack::kTPCrefit &&
696 track->GetStatus()&AliESDtrack::kITSrefit &&
697 track->GetTPCNcls() >=70 &&
698 track->GetStatus()&AliESDtrack::kTOFpid &&
699 track->GetStatus()&AliESDtrack::kTOFout &&
700 track->GetStatus()&AliESDtrack::kTIME)) continue;
701 AliAODPid *pid = track->GetDetPid();
a2121012 702 if(!pid) {if (fDebug>1)cout<<"No AliAODPid found"<<endl; continue;}
4e61a020 703
704 Double_t times[5];
705 pid->GetIntegratedTimes(times);
706
707 ((TH1F*)fChecks->FindObject("hTOFtime"))->Fill(times[3]);
708 ((TH1F*)fChecks->FindObject("hTOFtimeKaonHyptime"))->Fill(track->P(),pid->GetTOFsignal()-times[3]); //3 is kaon
709
710 ((TH1F*)fChecks->FindObject("hTOFsig"))->Fill(pid->GetTOFsignal());
6b3e3c78 711 ((TH1F*)fChecks->FindObject("hTPCsig"))->Fill(pid->GetTPCsignal());
712 ((TH1F*)fChecks->FindObject("hTPCsigvsp"))->Fill(track->P(),pid->GetTPCsignal());
4e61a020 713 if (pid->GetTOFsignal()< 0) ((TH1F*)fChecks->FindObject("hTOFsig"))->Fill(-1);
714
75638da0 715 //check clusters of the tracks
716 Int_t nclsTot=0,nclsSPD=0;
40445ada 717
75638da0 718 for(Int_t l=0;l<6;l++) {
719 if(TESTBIT(track->GetITSClusterMap(),l)) {
720 nclsTot++; if(l<2) nclsSPD++;
b272aebf 721 }
75638da0 722 }
6b3e3c78 723
75638da0 724 if (track->Pt()>0.3 &&
725 track->GetStatus()&AliESDtrack::kTPCrefit &&
726 track->GetStatus()&AliESDtrack::kITSrefit &&
727 nclsTot>3 &&
728 nclsSPD>0) {//fill hist good tracks
a2121012 729
75638da0 730 ((TH1F*)fChecks->FindObject("hptGoodTr"))->Fill(track->Pt());
731 isGoodTrack++;
40445ada 732 }
75638da0 733 //cout<<"isGoodTrack = "<<isGoodTrack<<endl;
734 ((TH1F*)fChecks->FindObject("hdistrGoodTr"))->Fill(isGoodTrack);
735 }
736 //number of events with good vertex and at least 2 good tracks
737 if (isGoodTrack>=2 && isGoodVtx) fNentries->Fill(4);
a96083b9 738
75638da0 739 // loop over candidates
740 Int_t nInD0toKpi = inputArray->GetEntriesFast();
4e61a020 741 if(fDebug>2) printf("Number of D0->Kpi: %d\n",nInD0toKpi);
75638da0 742
6b3e3c78 743 // FILE *f=fopen("4display.txt","a");
744 // fprintf(f,"Number of D0->Kpi: %d\n",nInD0toKpi);
a96083b9 745 Int_t nSelectedloose=0,nSelectedtight=0;
75638da0 746 for (Int_t iD0toKpi = 0; iD0toKpi < nInD0toKpi; iD0toKpi++) {
747 //Int_t nPosPairs=0, nNegPairs=0;
748 //cout<<"inside the loop"<<endl;
749 AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)inputArray->UncheckedAt(iD0toKpi);
4e61a020 750
751 //check daughters
752 if(!(d->GetDaughter(0) || d->GetDaughter(1))) {
753 AliDebug(1,"at least one daughter not found!");
754 fNentries->Fill(6);
a96083b9 755 continue;
4e61a020 756 }
757
6b3e3c78 758 // Bool_t unsetvtx=kFALSE;
759 // if(!d->GetOwnPrimaryVtx()) {
760 // d->SetOwnPrimaryVtx(vtx1); // needed to compute all variables
761 // unsetvtx=kTRUE;
762 // }
ea0d8716 763
4e61a020 764
75638da0 765 //check reco daughter in acceptance
4e61a020 766 /*
75638da0 767 Double_t eta0=d->EtaProng(0);
768 Double_t eta1=d->EtaProng(1);
4e61a020 769 */
770 if ( fCuts->IsInFiducialAcceptance(d->Pt(),d->Y(421)) ) {
771 //if( TMath::Abs(eta0)<0.9 && TMath::Abs(eta1)<0.9 ){
772 //apply cuts on tracks
a96083b9 773 nSelectedloose++;
3cc4604b 774 Int_t isSelected = fCuts->IsSelected(d,AliRDHFCuts::kTracks);
a2121012 775 if(((AliAODTrack*)d->GetDaughter(0))->GetTPCNcls() < 70 || ((AliAODTrack*)d->GetDaughter(1))->GetTPCNcls() < 70) isSelected=kFALSE;
a96083b9 776 if (!isSelected) continue;
777 nSelectedtight++;
4e61a020 778 fNentries->Fill(7);
779 if(fDebug>2) cout<<"tracks selected"<<endl;
3cc4604b 780
6b3e3c78 781 Int_t ptbin=fCuts->PtBin(d->Pt());
a96083b9 782 if(ptbin==-1) {fNentries->Fill(5); continue;} //out of bounds
4e61a020 783 FillVarHists(aod,d,mcArray,fCuts,fDistr);
6b3e3c78 784 FillMassHists(aod,d,mcArray,fCuts,fOutputMass);
75638da0 785 }
ea0d8716 786
6b3e3c78 787 //if(unsetvtx) d->UnsetOwnPrimaryVtx();
75638da0 788 } //end for prongs
a96083b9 789 fCounter->StoreCandidates(aod,nSelectedloose,kTRUE);
790 fCounter->StoreCandidates(aod,nSelectedtight,kFALSE);
40445ada 791 // Post the data
ea0d8716 792 PostData(1,fOutputMass);
793 PostData(2,fDistr);
40445ada 794 PostData(3,fNentries);
ea0d8716 795 PostData(4,fChecks);
a96083b9 796 PostData(6,fCounter);
40445ada 797 return;
798}
b272aebf 799
40445ada 800//____________________________________________________________________________
4e61a020 801void AliAnalysisTaskSED0Mass::FillVarHists(AliAODEvent* aod,AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC, AliRDHFCutsD0toKpi *cuts, TList *listout){
40445ada 802 //
803 // function used in UserExec to fill variable histograms:
804 //
b272aebf 805
40445ada 806 //add distr here
807 UInt_t pdgs[2];
808
809 Double_t mPDG=TDatabasePDG::Instance()->GetParticle(421)->Mass();
810 pdgs[0]=211;
811 pdgs[1]=321;
812 Double_t minvD0 = part->InvMassD0();
813 pdgs[1]=211;
814 pdgs[0]=321;
815 Double_t minvD0bar = part->InvMassD0bar();
816 //cout<<"inside mass cut"<<endl;
b272aebf 817
40445ada 818 Int_t pdgDgD0toKpi[2]={321,211};
819 Int_t lab=-9999;
820 if(fReadMC) lab=part->MatchToMC(421,arrMC,2,pdgDgD0toKpi); //return MC particle label if the array corresponds to a D0, -1 if not (cf. AliAODRecoDecay.cxx)
821 //Double_t pt = d->Pt(); //mother pt
6b3e3c78 822 Int_t isSelected=3;
823 Int_t isSelectedPID=99;
b272aebf 824
6b3e3c78 825 //isSelectedPID = cuts->IsSelected(part,AliRDHFCuts::kPID); //0 rejected,1 D0,2 Dobar, 3 both
826 isSelectedPID = cuts->IsSelectedPID(part); //0 rejected,1 D0,2 Dobar, 3 both
827 if (isSelectedPID==0)fNentries->Fill(8);
828 if (isSelectedPID==1)fNentries->Fill(9);
829 if (isSelectedPID==2)fNentries->Fill(10);
830 if (isSelectedPID==3)fNentries->Fill(11);
831 //fNentries->Fill(8+isSelectedPID);
3cc4604b 832
ea0d8716 833 if(fCutOnDistr){
6b3e3c78 834 isSelected = cuts->IsSelected(part,AliRDHFCuts::kCandidate,aod); //cuts with variables recalculated with new vertex (w/o daughters)
ea0d8716 835 if (!isSelected){
836 //cout<<"Not Selected"<<endl;
837 return;
838 }
839 }
0108fa62 840
4e61a020 841 Double_t invmasscut=0.03;
842
40445ada 843 TString fillthispi="",fillthisK="",fillthis="";
b272aebf 844
ea0d8716 845 Int_t ptbin=cuts->PtBin(part->Pt());
6b3e3c78 846
4e61a020 847
848 //recalculate vertex
849 AliAODVertex *vtxProp=0x0;
850 vtxProp=GetPrimaryVtxSkipped(aod,part);
851 Double_t dz1[2],dz2[2],covar1[3],covar2[3];//,d0xd0proper,errd0xd0proper;
852 dz1[0]=-99; dz2[0]=-99;
853 if(vtxProp) {
854 part->SetOwnPrimaryVtx(vtxProp);
855 //Bool_t unsetvtx=kTRUE;
856 //Calculate d0 for daughter tracks with Vtx Skipped
857 AliESDtrack *esdtrack1=new AliESDtrack((AliVTrack*)part->GetDaughter(0));
858 AliESDtrack *esdtrack2=new AliESDtrack((AliVTrack*)part->GetDaughter(1));
859 esdtrack1->PropagateToDCA(vtxProp,aod->GetMagneticField(),1.,dz1,covar1);
860 esdtrack2->PropagateToDCA(vtxProp,aod->GetMagneticField(),1.,dz2,covar2);
6b3e3c78 861 delete vtxProp; vtxProp=NULL;
4e61a020 862 delete esdtrack1;
863 delete esdtrack2;
864 }
865
866 Double_t d0[2]={dz1[0],dz2[0]};
6b3e3c78 867 Double_t decl[2]={part->DecayLength(),part->NormalizedDecayLength()};
868 part->UnsetOwnPrimaryVtx();
4e61a020 869
40445ada 870 if(!fCutOnDistr || (fCutOnDistr && isSelected)){ //if no cuts or cuts passed
ea0d8716 871 //printf("\nif no cuts or cuts passed\n");
0108fa62 872
6b3e3c78 873 //disable the PID
874 if(!fUsePid4Distr) isSelectedPID=0;
875 if((lab>=0 && fReadMC) || (!fReadMC && isSelectedPID)){ //signal (from MC or PID)
876 //if(!vtxProp) fNentries->Fill(3);
40445ada 877 //check pdg of the prongs
878 AliAODTrack *prong0=(AliAODTrack*)part->GetDaughter(0);
879 AliAODTrack *prong1=(AliAODTrack*)part->GetDaughter(1);
4e61a020 880 if(!prong0 || !prong1) {
881 return;
882 }
40445ada 883 Int_t labprong[2];
6b3e3c78 884 if(fReadMC){
885 labprong[0]=prong0->GetLabel();
886 labprong[1]=prong1->GetLabel();
887 }
40445ada 888 AliAODMCParticle *mcprong=0;
889 Int_t pdgProng[2]={0,0};
6b3e3c78 890
40445ada 891 for (Int_t iprong=0;iprong<2;iprong++){
6b3e3c78 892 if(fReadMC && labprong[iprong]>=0) {
893 mcprong= (AliAODMCParticle*)arrMC->At(labprong[iprong]);
894 pdgProng[iprong]=mcprong->GetPdgCode();
895 }
0108fa62 896 }
6b3e3c78 897
0108fa62 898
40445ada 899 //no mass cut ditributions: ptbis
900
901 fillthispi="hptpiSnoMcut_";
902 fillthispi+=ptbin;
903
904 fillthisK="hptKSnoMcut_";
905 fillthisK+=ptbin;
906
6b3e3c78 907 if ((TMath::Abs(pdgProng[0]) == 211 && TMath::Abs(pdgProng[1]) == 321)
908 || (isSelectedPID==1 || isSelectedPID==3)){
40445ada 909 ((TH1F*)listout->FindObject(fillthispi))->Fill(part->PtProng(0));
910 ((TH1F*)listout->FindObject(fillthisK))->Fill(part->PtProng(1));
6b3e3c78 911 }
912
913 if ((TMath::Abs(pdgProng[0]) == 321 && TMath::Abs(pdgProng[1]) == 211)
914 || (isSelectedPID==2 || isSelectedPID==3)){
915 ((TH1F*)listout->FindObject(fillthisK))->Fill(part->PtProng(0));
916 ((TH1F*)listout->FindObject(fillthispi))->Fill(part->PtProng(1));
40445ada 917 }
918
919 //no mass cut ditributions: mass
920 fillthis="hMassS_";
921 fillthis+=ptbin;
922
6b3e3c78 923 if ((fReadMC && ((AliAODMCParticle*)arrMC->At(lab))->GetPdgCode() == 421)
924 || (!fReadMC && (isSelectedPID==1 || isSelectedPID==3))){//D0
40445ada 925 ((TH1F*)listout->FindObject(fillthis))->Fill(minvD0);
926 }
927 else { //D0bar
6b3e3c78 928 if(fReadMC || (!fReadMC && isSelectedPID > 1))
40445ada 929 ((TH1F*)listout->FindObject(fillthis))->Fill(minvD0bar);
930 }
4e61a020 931
0108fa62 932 //apply cut on invariant mass on the pair
4e61a020 933 if(TMath::Abs(minvD0-mPDG)<invmasscut || TMath::Abs(minvD0bar-mPDG)<invmasscut){
40445ada 934
0108fa62 935 if(fArray==1) cout<<"LS signal: ERROR"<<endl;
936 for (Int_t iprong=0; iprong<2; iprong++){
40445ada 937 AliAODTrack *prong=(AliAODTrack*)part->GetDaughter(iprong);
6b3e3c78 938 if (fReadMC) labprong[iprong]=prong->GetLabel();
40445ada 939
0108fa62 940 //cout<<"prong name = "<<prong->GetName()<<" label = "<<prong->GetLabel()<<endl;
6b3e3c78 941 Int_t pdgprong=0;
942 if(fReadMC && labprong[iprong]>=0) {
943 mcprong= (AliAODMCParticle*)arrMC->At(labprong[iprong]);
944 pdgprong=mcprong->GetPdgCode();
945 }
946
947 Bool_t isPionHere[2]={(isSelectedPID==1 || isSelectedPID==3) ? kTRUE : kFALSE,(isSelectedPID==2 || isSelectedPID==3) ? kTRUE : kFALSE};
948
949 if(TMath::Abs(pdgprong)==211 || isPionHere[iprong]) {
0108fa62 950 //cout<<"pi"<<endl;
40445ada 951
952 fillthispi="hptpiS_";
953 fillthispi+=ptbin;
954 ((TH1F*)listout->FindObject(fillthispi))->Fill(part->PtProng(iprong));
4e61a020 955 fillthispi="hd0piS_";
956 fillthispi+=ptbin;
957 ((TH1F*)listout->FindObject(fillthispi))->Fill(part->Getd0Prong(iprong));
958 if(d0[iprong] != -99) {
959
960 fillthispi="hd0vpiS_";
961 fillthispi+=ptbin;
962 ((TH1F*)listout->FindObject(fillthispi))->Fill(d0[iprong]);
963 }
964
b272aebf 965 }
40445ada 966
6b3e3c78 967 if(TMath::Abs(pdgprong)==321 || !isPionHere[iprong]) {
40445ada 968 //cout<<"kappa"<<endl;
969
970 fillthisK="hptKS_";
971 fillthisK+=ptbin;
972 ((TH1F*)listout->FindObject(fillthisK))->Fill(part->PtProng(iprong));
40445ada 973 fillthisK="hd0KS_";
974 fillthisK+=ptbin;
975 ((TH1F*)listout->FindObject(fillthisK))->Fill(part->Getd0Prong(iprong));
4e61a020 976 if (d0[iprong] != -99){
977 fillthisK="hd0vKS_";
978 fillthisK+=ptbin;
979 ((TH1F*)listout->FindObject(fillthisK))->Fill(d0[iprong]);
980 }
b272aebf 981 }
b272aebf 982
4e61a020 983 fillthis="hcosthpointd0S_";
40445ada 984 fillthis+=ptbin;
4e61a020 985 ((TH1F*)listout->FindObject(fillthis))->Fill(part->CosPointingAngle(),part->Getd0Prong(iprong));
b272aebf 986
4e61a020 987 } //end loop on prongs
b272aebf 988
4e61a020 989 fillthis="hdcaS_";
990 fillthis+=ptbin;
991 ((TH1F*)listout->FindObject(fillthis))->Fill(part->GetDCA());
992
993 fillthis="hcosthetapointS_";
994 fillthis+=ptbin;
995 ((TH1F*)listout->FindObject(fillthis))->Fill(part->CosPointingAngle());
996
997 fillthis="hcosthpointd0d0S_";
998 fillthis+=ptbin;
999 ((TH2F*)listout->FindObject(fillthis))->Fill(part->CosPointingAngle(),part->Prodd0d0());
b272aebf 1000
4e61a020 1001 fillthis="hcosthetastarS_";
1002 fillthis+=ptbin;
6b3e3c78 1003 if ((fReadMC && ((AliAODMCParticle*)arrMC->At(lab))->GetPdgCode() == 421)) ((TH1F*)listout->FindObject(fillthis))->Fill(part->CosThetaStarD0());
1004 else {
1005 if (fReadMC || isSelectedPID>1)((TH1F*)listout->FindObject(fillthis))->Fill(part->CosThetaStarD0bar());
1006 if(isSelectedPID==1 || isSelectedPID==3)((TH1F*)listout->FindObject(fillthis))->Fill(part->CosThetaStarD0());
1007 }
4e61a020 1008 fillthis="hd0d0S_";
1009 fillthis+=ptbin;
1010 ((TH1F*)listout->FindObject(fillthis))->Fill(part->Prodd0d0());
1011
1012 if(d0[0] != -99 && d0[1] != -99 ){
1013 fillthis="hd0d0vS_";
40445ada 1014 fillthis+=ptbin;
4e61a020 1015 ((TH1F*)listout->FindObject(fillthis))->Fill(d0[0]*d0[1]);
b272aebf 1016 }
4e61a020 1017
1018 fillthis="hdeclS_";
1019 fillthis+=ptbin;
1020 ((TH1F*)listout->FindObject(fillthis))->Fill(part->DecayLength());
1021
1022 fillthis="hnormdeclS_";
1023 fillthis+=ptbin;
1024 ((TH1F*)listout->FindObject(fillthis))->Fill(part->NormalizedDecayLength());
a2121012 1025
1026 if(vtxProp) {
1027 fillthis="hdeclvS_";
1028 fillthis+=ptbin;
6b3e3c78 1029 ((TH1F*)listout->FindObject(fillthis))->Fill(decl[0]);
a2121012 1030
1031 fillthis="hnormdeclvS_";
1032 fillthis+=ptbin;
6b3e3c78 1033 ((TH1F*)listout->FindObject(fillthis))->Fill(decl[1]);
a2121012 1034 }
6b3e3c78 1035 } //end mass cut
4e61a020 1036
0108fa62 1037 } else{ //Background or LS
6b3e3c78 1038 //if(!fReadMC){
0108fa62 1039 //cout<<"is background"<<endl;
40445ada 1040
b272aebf 1041 //no mass cut distributions: mass, ptbis
40445ada 1042 fillthis="hMassB_";
1043 fillthis+=ptbin;
ea0d8716 1044
1045 if (!fCutOnDistr || (fCutOnDistr && (isSelected==1 || isSelected==3))) ((TH1F*)listout->FindObject(fillthis))->Fill(minvD0);
1046 if (!fCutOnDistr || (fCutOnDistr && isSelected>1)) ((TH1F*)listout->FindObject(fillthis))->Fill(minvD0bar);
b272aebf 1047
40445ada 1048 fillthis="hptB1prongnoMcut_";
1049 fillthis+=ptbin;
1050
1051 ((TH1F*)listout->FindObject(fillthis))->Fill(part->PtProng(0));
1052
1053 fillthis="hptB2prongsnoMcut_";
1054 fillthis+=ptbin;
1055 ((TH1F*)listout->FindObject(fillthis))->Fill(part->PtProng(0));
1056 ((TH1F*)listout->FindObject(fillthis))->Fill(part->PtProng(1));
0108fa62 1057
1058 //apply cut on invariant mass on the pair
4e61a020 1059 if(TMath::Abs(minvD0-mPDG)<invmasscut || TMath::Abs(minvD0bar-mPDG)<invmasscut){
0108fa62 1060
1061
6b3e3c78 1062 AliAODTrack *prongg=(AliAODTrack*)part->GetDaughter(0);
1063 if(!prongg) {
4e61a020 1064 if(fDebug>2) cout<<"No daughter found";
1065 return;
1066 }
0108fa62 1067 else{
449b1302 1068 if(prongg->Charge()==1) {
1069 fTotPosPairs[ptbin]++;
1070 ((TH1F*)fDistr->FindObject("hpospair"))->Fill(ptbin);
1071 } else {
1072 fTotNegPairs[ptbin]++;
1073 ((TH1F*)fDistr->FindObject("hnegpair"))->Fill(ptbin);
1074 }
0108fa62 1075 }
40445ada 1076
1077 //normalise pt distr to half afterwards
1078 fillthis="hptB_";
1079 fillthis+=ptbin;
4e61a020 1080 ((TH1F*)listout->FindObject(fillthis))->Fill(part->PtProng(0));
1081 ((TH1F*)listout->FindObject(fillthis))->Fill(part->PtProng(1));
0108fa62 1082
4e61a020 1083 fillthis="hd0p0B_";
40445ada 1084 fillthis+=ptbin;
1085 ((TH1F*)listout->FindObject(fillthis))->Fill(part->Getd0Prong(0));
4e61a020 1086 fillthis="hd0p1B_";
1087 fillthis+=ptbin;
1088 ((TH1F*)listout->FindObject(fillthis))->Fill(part->Getd0Prong(1));
a2121012 1089 fillthis="hd0B_";
1090 fillthis+=ptbin;
1091 ((TH1F*)listout->FindObject(fillthis))->Fill(part->Getd0Prong(0));
1092 ((TH1F*)listout->FindObject(fillthis))->Fill(part->Getd0Prong(1));
4e61a020 1093
6b3e3c78 1094 if(fReadMC){
1095 Int_t pdgMother[2]={0,0};
1096 Double_t factor[2]={1,1};
1097
1098 for(Int_t iprong=0;iprong<2;iprong++){
1099 AliAODTrack *prong=(AliAODTrack*)part->GetDaughter(iprong);
1100 lab=prong->GetLabel();
1101 if(lab>=0){
1102 AliAODMCParticle* mcprong=(AliAODMCParticle*)arrMC->At(lab);
1103 if(mcprong){
1104 Int_t labmom=mcprong->GetMother();
1105 if(labmom>=0){
1106 AliAODMCParticle* mcmother=(AliAODMCParticle*)arrMC->At(labmom);
1107 if(mcmother) pdgMother[iprong]=mcmother->GetPdgCode();
1108 }
1109 }
1110 }
1111
1112 fillthis="hd0moresB_";
1113 fillthis+=ptbin;
1114
1115 if(TMath::Abs(pdgMother[iprong])==310 || TMath::Abs(pdgMother[iprong])==130 || TMath::Abs(pdgMother[iprong])==321){ //K^0_S, K^0_L, K^+-
1116 if(part->PtProng(iprong)<=1)factor[iprong]=1./.7;
1117 else factor[iprong]=1./.6;
1118 fNentries->Fill(12);
1119 }
1120
1121 if(TMath::Abs(pdgMother[iprong])==3122) { //Lambda
1122 factor[iprong]=1./0.25;
1123 fNentries->Fill(13);
1124 }
1125 fillthis="hd0moresB_";
1126 fillthis+=ptbin;
1127
1128 ((TH1F*)listout->FindObject(fillthis))->Fill(part->Getd0Prong(iprong),factor[iprong]);
1129
1130 fillthis="hd0vmoresB_";
1131 fillthis+=ptbin;
1132 ((TH1F*)listout->FindObject(fillthis))->Fill(d0[iprong],factor[iprong]);
1133
1134 }
1135
1136 fillthis="hd0d0moresB_";
1137 fillthis+=ptbin;
1138 ((TH1F*)listout->FindObject(fillthis))->Fill(part->Prodd0d0(),factor[0]*factor[1]);
1139
1140 fillthis="hd0d0vmoresB_";
1141 fillthis+=ptbin;
1142 ((TH1F*)listout->FindObject(fillthis))->Fill(d0[0]*d0[1],factor[0]*factor[1]);
1143
1144 fillthis="hcosthetapointmoresB_";
1145 fillthis+=ptbin;
1146 ((TH1F*)listout->FindObject(fillthis))->Fill(part->CosPointingAngle(),factor[0]*factor[1]);
1147 }
1148
4e61a020 1149 fillthis="hd0vp0B_";
1150 fillthis+=ptbin;
1151 ((TH1F*)listout->FindObject(fillthis))->Fill(d0[0]);
1152 fillthis="hd0vp1B_";
1153 fillthis+=ptbin;
1154 ((TH1F*)listout->FindObject(fillthis))->Fill(d0[1]);
1155
a2121012 1156 fillthis="hd0vB_";
1157 fillthis+=ptbin;
1158 ((TH1F*)listout->FindObject(fillthis))->Fill(d0[0]);
1159 ((TH1F*)listout->FindObject(fillthis))->Fill(d0[1]);
1160
1161
4e61a020 1162 fillthis="hcosthpointd0B_";
1163 fillthis+=ptbin;
1164 ((TH1F*)listout->FindObject(fillthis))->Fill(part->CosPointingAngle(),part->Getd0Prong(0));
1165 ((TH1F*)listout->FindObject(fillthis))->Fill(part->CosPointingAngle(),part->Getd0Prong(1));
1166
ce39f0ac 1167
40445ada 1168 fillthis="hdcaB_";
1169 fillthis+=ptbin;
1170 ((TH1F*)listout->FindObject(fillthis))->Fill(part->GetDCA());
34137226 1171
40445ada 1172 fillthis="hcosthetastarB_";
1173 fillthis+=ptbin;
ea0d8716 1174 if (!fCutOnDistr || (fCutOnDistr && (isSelected==1 || isSelected==3)))((TH1F*)listout->FindObject(fillthis))->Fill(part->CosThetaStarD0());
4e61a020 1175 if (!fCutOnDistr || (fCutOnDistr && isSelected>1))((TH1F*)listout->FindObject(fillthis))->Fill(part->CosThetaStarD0bar());
40445ada 1176
1177 fillthis="hd0d0B_";
1178 fillthis+=ptbin;
1179 ((TH1F*)listout->FindObject(fillthis))->Fill(part->Prodd0d0());
1180
a2121012 1181 if(d0[0] != -99 && d0[1] != -99 ){
1182 fillthis="hd0d0vB_";
1183 fillthis+=ptbin;
1184 ((TH1F*)listout->FindObject(fillthis))->Fill(d0[0]*d0[1]);
1185 }
1186
40445ada 1187 fillthis="hcosthetapointB_";
1188 fillthis+=ptbin;
1189 ((TH1F*)listout->FindObject(fillthis))->Fill(part->CosPointingAngle());
1190
1191 fillthis="hcosthpointd0d0B_";
1192 fillthis+=ptbin;
1193 ((TH2F*)listout->FindObject(fillthis))->Fill(part->CosPointingAngle(),part->Prodd0d0());
a2121012 1194
1195 fillthis="hdeclB_";
1196 fillthis+=ptbin;
1197 ((TH1F*)listout->FindObject(fillthis))->Fill(part->DecayLength());
1198
1199 fillthis="hnormdeclB_";
1200 fillthis+=ptbin;
1201 ((TH1F*)listout->FindObject(fillthis))->Fill(part->NormalizedDecayLength());
1202
1203 if(vtxProp) {
1204 fillthis="hdeclvB_";
1205 fillthis+=ptbin;
6b3e3c78 1206 ((TH1F*)listout->FindObject(fillthis))->Fill(decl[0]);
a2121012 1207
1208 fillthis="hnormdeclvB_";
1209 fillthis+=ptbin;
6b3e3c78 1210 ((TH1F*)listout->FindObject(fillthis))->Fill(decl[1]);
a2121012 1211 }
6b3e3c78 1212 }//mass cut
4e61a020 1213 }//else (background)
40445ada 1214 }
49061176 1215 return;
1216}
ea0d8716 1217//____________________________________________________________________________
1218
6b3e3c78 1219void AliAnalysisTaskSED0Mass::FillMassHists(AliAODEvent* aod,AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC, AliRDHFCutsD0toKpi* cuts, TList *listout){
49061176 1220 //
40445ada 1221 // function used in UserExec to fill mass histograms:
49061176 1222 //
9de8c723 1223
1224
1225 Double_t mPDG=TDatabasePDG::Instance()->GetParticle(421)->Mass();
1226
6b3e3c78 1227 Int_t isSelected=cuts->IsSelected(part,AliRDHFCuts::kCandidate,aod); //selected
1228
ea0d8716 1229 //cout<<"is selected = "<<isSelected<<endl;
feb73eca 1230
ea0d8716 1231 //cout<<"check cuts = "<<endl;
1232 //cuts->PrintAll();
1233 if (!isSelected){
1234 //cout<<"Not Selected"<<endl;
6b3e3c78 1235 //cout<<"Rejected because "<<cuts->GetWhy()<<endl;
ea0d8716 1236 return;
1237 }
9de8c723 1238
6b3e3c78 1239
4e61a020 1240 if(fDebug>2) cout<<"Candidate selected"<<endl;
a41f6fad 1241
ea0d8716 1242 Double_t invmassD0 = part->InvMassD0(), invmassD0bar = part->InvMassD0bar();
1243 //printf("SELECTED\n");
1244 Int_t ptbin=cuts->PtBin(part->Pt());
9de8c723 1245
ea0d8716 1246 AliAODTrack *prong=(AliAODTrack*)part->GetDaughter(0);
4e61a020 1247 if(!prong) {
1248 AliDebug(2,"No daughter found");
1249 return;
1250 }
ea0d8716 1251 else{
4b4a1d25 1252 if(prong->Charge()==1) {fTotPosPairs[ptbin]++;} else {fTotNegPairs[ptbin]++;}
ea0d8716 1253 }
7646d6da 1254
6b3e3c78 1255 for(Int_t it=0;it<2;it++){
1256 //fill TOFsig-kaontime
1257 AliAODTrack* track=(AliAODTrack*)part->GetDaughter(it);
1258 AliAODPid *pid = track->GetDetPid();
1259 Double_t times[5];
1260 pid->GetIntegratedTimes(times);
1261 ((TH2F*)fChecks->FindObject("hTOFtimeKaonHyptimeAC"))->Fill(track->P(),pid->GetTOFsignal()-times[3]);
1262 }
1263
1264 //request on spd points to be addes
1265 if(/*nSPD==2 && */part->Pt() > 5. && (TMath::Abs(invmassD0-mPDG)<0.01 || TMath::Abs(invmassD0bar-mPDG)<0.01)){
1266 FILE *f=fopen("4display.txt","a");
1267 fprintf(f,"pt: %f \n Rapidity: %f \t Period Number: %x \t Run Number: %d \t BunchCrossNumb: %d \t OrbitNumber: %d\n",part->Pt(),part->Y(421),aod->GetPeriodNumber(),aod->GetRunNumber(),aod->GetBunchCrossNumber(),aod->GetOrbitNumber());
1268 fclose(f);
1269 //printf("PrimVtx NContributors: %d \n Prongs Rel Angle: %f \n \n",ncont,relangle);
1270 }
1271
1272
1273
ea0d8716 1274 TString fillthis="";
1275 Int_t pdgDgD0toKpi[2]={321,211};
1276 Int_t labD0=-1;
1277 if (fReadMC) labD0 = part->MatchToMC(421,arrMC,2,pdgDgD0toKpi); //return MC particle label if the array corresponds to a D0, -1 if not (cf. AliAODRecoDecay.cxx)
1278
1279 //count candidates selected by cuts
1280 fNentries->Fill(1);
1281 //count true D0 selected by cuts
1282 if (fReadMC && labD0>=0) fNentries->Fill(2);
1283 //PostData(3,fNentries);
1284
449b1302 1285 if ((isSelected==1 || isSelected==3) && fFillOnlyD0D0bar<2) { //D0
ea0d8716 1286 fillthis="histMass_";
1287 fillthis+=ptbin;
1288 //cout<<"Filling "<<fillthis<<endl;
1289
1290 //printf("Fill mass with D0");
1291 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
feb73eca 1292
ea0d8716 1293 if(labD0>=0) {
1294 if(fArray==1) cout<<"LS signal ERROR"<<endl;
1295
1296 AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
1297 Int_t pdgD0 = partD0->GetPdgCode();
1298 //cout<<"pdg = "<<pdgD0<<endl;
1299 if (pdgD0==421){ //D0
1300 //cout<<"Fill S with D0"<<endl;
1301 fillthis="histSgn_";
1302 fillthis+=ptbin;
1303 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
1304 if(TMath::Abs(invmassD0 - mPDG) < 0.027){
1305 fillthis="histSgn27_";
46c96ce5 1306 fillthis+=ptbin;
1307 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
a4ae02cd 1308 }
ea0d8716 1309 } else{ //it was a D0bar
1310 fillthis="histRfl_";
a4ae02cd 1311 fillthis+=ptbin;
1312 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
1313 }
ea0d8716 1314 } else {//background
1315 fillthis="histBkg_";
a4ae02cd 1316 fillthis+=ptbin;
ea0d8716 1317 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0);
1318 }
1319
1320 }
449b1302 1321 if (isSelected>1 && (fFillOnlyD0D0bar==0 || fFillOnlyD0D0bar==2)) { //D0bar
ea0d8716 1322 fillthis="histMass_";
1323 fillthis+=ptbin;
1324 //printf("Fill mass with D0bar");
1325 ((TH1F*)listout->FindObject(fillthis))->Fill(invmassD0bar);
feb73eca 1326
ea0d8716 1327 if(labD0>=0) {
1328 if(fArray==1) cout<<"LS signal ERROR"<<endl;
1329 AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
1330 Int_t pdgD0 = partD0->GetPdgCode();
1331 //cout<<" pdg = "<<pdgD0<<endl;
1332 if (pdgD0==-421){ //D0bar
1333 fillthis="histSgn_";
1334 fillthis+=ptbin;
1335 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
1336 if (TMath::Abs(invmassD0bar - mPDG) < 0.027){
1337 fillthis="histSgn27_";
a4ae02cd 1338 fillthis+=ptbin;
1339 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
ea0d8716 1340 }
9de8c723 1341
a4ae02cd 1342
ea0d8716 1343 } else{
1344 fillthis="histRfl_";
a4ae02cd 1345 fillthis+=ptbin;
1346 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
1347 }
ea0d8716 1348 } else {//background or LS
1349 fillthis="histBkg_";
1350 fillthis+=ptbin;
1351 ((TH1F*)(listout->FindObject(fillthis)))->Fill(invmassD0bar);
a4ae02cd 1352 }
ea0d8716 1353 }
a4ae02cd 1354
40445ada 1355 return;
49061176 1356}
4e61a020 1357
1358//__________________________________________________________________________
1359AliAODVertex* AliAnalysisTaskSED0Mass::GetPrimaryVtxSkipped(AliAODEvent *aodev,AliAODRecoDecayHF2Prong *d){
1360 //Calculate the primary vertex w/o the daughter tracks of the candidate
1361
4e61a020 1362 Int_t skipped[2];
1363 Int_t nTrksToSkip=2;
1364 AliAODTrack *dgTrack = (AliAODTrack*)d->GetDaughter(0);
1365 if(!dgTrack){
1366 AliDebug(2,"no daughter found!");
1367 return 0x0;
1368 }
1369 skipped[0]=dgTrack->GetID();
1370 dgTrack = (AliAODTrack*)d->GetDaughter(1);
1371 if(!dgTrack){
1372 AliDebug(2,"no daughter found!");
1373 return 0x0;
1374 }
1375 skipped[1]=dgTrack->GetID();
1376
6b3e3c78 1377 AliESDVertex *vertexESD=0x0;
1378 AliAODVertex *vertexAOD=0x0;
1379 AliVertexerTracks *vertexer = new AliVertexerTracks(aodev->GetMagneticField());
1380
4e61a020 1381 //
1382 vertexer->SetSkipTracks(nTrksToSkip,skipped);
4e61a020 1383 vertexer->SetMinClusters(4);
a2121012 1384 vertexESD = (AliESDVertex*)vertexer->FindPrimaryVertex(aodev);
4e61a020 1385 if(!vertexESD) return vertexAOD;
1386 if(vertexESD->GetNContributors()<=0) {
1387 AliDebug(2,"vertexing failed");
1388 delete vertexESD; vertexESD=NULL;
1389 return vertexAOD;
1390 }
1391
1392 delete vertexer; vertexer=NULL;
1393
1394
1395 // convert to AliAODVertex
1396 Double_t pos[3],cov[6],chi2perNDF;
1397 vertexESD->GetXYZ(pos); // position
1398 vertexESD->GetCovMatrix(cov); //covariance matrix
1399 chi2perNDF = vertexESD->GetChi2toNDF();
1400 delete vertexESD; vertexESD=NULL;
1401
1402 vertexAOD = new AliAODVertex(pos,cov,chi2perNDF);
1403 return vertexAOD;
1404
1405}
1406
1407
49061176 1408//________________________________________________________________________
1409void AliAnalysisTaskSED0Mass::Terminate(Option_t */*option*/)
1410{
1411 // Terminate analysis
1412 //
1413 if(fDebug > 1) printf("AnalysisTaskSED0Mass: Terminate() \n");
1414
6321ee46 1415
ea0d8716 1416 fOutputMass = dynamic_cast<TList*> (GetOutputData(1));
1417 if (!fOutputMass) {
1418 printf("ERROR: fOutputMass not available\n");
a4ae02cd 1419 return;
1420 }
ea0d8716 1421 fDistr = dynamic_cast<TList*> (GetOutputData(2));
1422 if (!fDistr) {
1423 printf("ERROR: fDistr not available\n");
a41f6fad 1424 return;
1425 }
5b2e5fae 1426
40445ada 1427 fNentries = dynamic_cast<TH1F*>(GetOutputData(3));
5b2e5fae 1428
40445ada 1429 if(!fNentries){
1430 printf("ERROR: fNEntries not available\n");
1431 return;
1432 }
5b2e5fae 1433
ea0d8716 1434 fChecks = dynamic_cast<TList*> (GetOutputData(4));
34137226 1435 if (!fChecks) {
1436 printf("ERROR: fChecks not available\n");
1437 return;
1438 }
a96083b9 1439 fCounter = dynamic_cast<AliNormalizationCounter*>(GetOutputData(6));
1440 if (!fCounter) {
1441 printf("ERROR: fCounter not available\n");
1442 return;
1443 }
4e61a020 1444 for(Int_t ipt=0;ipt<5;ipt++){ //change 5 in GetNPtBins when sure it is written and check
1445
4e61a020 1446
1447 if(fArray==1){
feb73eca 1448 fLsNormalization = 2.*TMath::Sqrt(fTotPosPairs[ipt]*fTotNegPairs[ipt]);
1449
1450
ea0d8716 1451 if(fLsNormalization>1e-6) {
9de8c723 1452
feb73eca 1453 TString massName="histMass_";
ea0d8716 1454 massName+=ipt;
1455 ((TH1F*)fOutputMass->FindObject(massName))->Scale((1/fLsNormalization)*((TH1F*)fOutputMass->FindObject(massName))->GetEntries());
1456
feb73eca 1457 }
40445ada 1458
feb73eca 1459
40445ada 1460 fLsNormalization = 2.*TMath::Sqrt(fTotPosPairs[4]*fTotNegPairs[4]);
1461
ea0d8716 1462 if(fLsNormalization>1e-6) {
40445ada 1463
1464 TString nameDistr="hptB_";
ea0d8716 1465 nameDistr+=ipt;
40445ada 1466 ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
1467 nameDistr="hdcaB_";
ea0d8716 1468 nameDistr+=ipt;
40445ada 1469 ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
1470 nameDistr="hcosthetastarB_";
ea0d8716 1471 nameDistr+=ipt;
40445ada 1472 ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
1473 nameDistr="hd0B_";
ea0d8716 1474 nameDistr+=ipt;
40445ada 1475 ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
1476 nameDistr="hd0d0B_";
ea0d8716 1477 nameDistr+=ipt;
40445ada 1478 ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
1479 nameDistr="hcosthetapointB_";
ea0d8716 1480 nameDistr+=ipt;
40445ada 1481 ((TH1F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH1F*)fDistr->FindObject(nameDistr))->GetEntries());
1482 nameDistr="hcosthpointd0d0B_";
ea0d8716 1483 nameDistr+=ipt;
40445ada 1484 ((TH2F*)fDistr->FindObject(nameDistr))->Scale((1/fLsNormalization)*((TH2F*)fDistr->FindObject(nameDistr))->GetEntries());
feb73eca 1485
40445ada 1486 }
feb73eca 1487 }
1488 }
4e61a020 1489 TString cvname,cstname;
527f330b 1490
1491 if (fArray==0){
1492 cvname="D0invmass";
4e61a020 1493 cstname="cstat0";
1494 } else {
1495 cvname="LSinvmass";
1496 cstname="cstat1";
1497 }
527f330b 1498
34137226 1499 TCanvas *cMass=new TCanvas(cvname,cvname);
1500 cMass->cd();
ea0d8716 1501 ((TH1F*)fOutputMass->FindObject("histMass_3"))->Draw();
527f330b 1502
4e61a020 1503 TCanvas* cStat=new TCanvas(cstname,Form("Stat%s",fArray ? "LS" : "D0"));
34137226 1504 cStat->cd();
1505 cStat->SetGridy();
1506 fNentries->Draw("htext0");
527f330b 1507
4e61a020 1508 // TCanvas *ccheck=new TCanvas(Form("cc%d",fArray),Form("cc%d",fArray));
1509 // ccheck->cd();
1510
49061176 1511 return;
1512}
1513