]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/vertexingHF/AliAnalysisTaskSELc2V0bachelorTMVA.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'master' of https://git.cern.ch/reps/AliRoot
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliAnalysisTaskSELc2V0bachelorTMVA.cxx
CommitLineData
0abe5247 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 * appeuear 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/* $Id: AliAnalysisTaskSELc2V0bachelorTMVA.cxx 62231 2013-04-29 09:47:26Z fprino $ */
17
18//
19//
20// Base class for Lc2V0 Analysis to be used with TMVA
21//
22//
23
24//------------------------------------------------------------------------------------------
25//
26// Author: C. Zampolli (from AliAnalysisTaskSELc2V0bachelor class by A.De Caro, P. Pagano)
27//
28//------------------------------------------------------------------------------------------
29
30#include <TSystem.h>
31#include <TParticle.h>
32#include <TParticlePDG.h>
33#include <TH1F.h>
34#include <TH1I.h>
35#include <TH2F.h>
36#include <TTree.h>
37#include "TROOT.h"
38#include <TDatabasePDG.h>
ef87f0d6 39#include "TVector3.h"
40
0abe5247 41#include <AliAnalysisDataSlot.h>
42#include <AliAnalysisDataContainer.h>
43#include "AliStack.h"
44#include "AliMCEvent.h"
45#include "AliAnalysisManager.h"
46#include "AliAODMCHeader.h"
47#include "AliAODHandler.h"
48#include "AliLog.h"
49#include "AliAODVertex.h"
50#include "AliAODRecoDecay.h"
51#include "AliAODRecoDecayHF.h"
52#include "AliAODRecoCascadeHF.h"
53#include "AliAnalysisVertexingHF.h"
54#include "AliESDtrack.h"
55#include "AliAODTrack.h"
56#include "AliAODv0.h"
57#include "AliAODMCParticle.h"
58#include "AliAnalysisTaskSE.h"
59#include "AliAnalysisTaskSELc2V0bachelorTMVA.h"
60#include "AliNormalizationCounter.h"
61#include "AliAODPidHF.h"
62#include "AliPIDResponse.h"
63#include "AliPIDCombined.h"
64#include "AliTOFPIDResponse.h"
65#include "AliInputEventHandler.h"
66#include "AliAODRecoDecayHF3Prong.h"
67#include "AliKFParticle.h"
68#include "AliKFVertex.h"
ef87f0d6 69#include "AliExternalTrackParam.h"
70#include "AliESDtrack.h"
0abe5247 71
72using std::cout;
73using std::endl;
74
75ClassImp(AliAnalysisTaskSELc2V0bachelorTMVA)
76
77//__________________________________________________________________________
78AliAnalysisTaskSELc2V0bachelorTMVA::AliAnalysisTaskSELc2V0bachelorTMVA():
79AliAnalysisTaskSE(),
80 fUseMCInfo(kFALSE),
81 fOutput(0),
82 fCEvents(0),
83 fPIDResponse(0),
84 fPIDCombined(0),
85 fIsK0sAnalysis(kFALSE),
86 fCounter(0),
87 fAnalCuts(0),
88 fListCuts(0),
89 fUseOnTheFlyV0(kFALSE),
90 fIsEventSelected(kFALSE),
91 fVariablesTreeSgn(0),
92 fVariablesTreeBkg(0),
93 fCandidateVariables(),
94 fIspA(kFALSE),
95 fHistoEvents(0),
96 fHistoLc(0),
97 fHistoLcOnTheFly(0),
98 fFillOnlySgn(kFALSE),
ef87f0d6 99 fHistoLcBeforeCuts(0),
100 fHistoFiducialAcceptance(0),
101 fHistoCodesSgn(0),
102 fHistoCodesBkg(0),
0abe5247 103 fHistoLcpKpiBeforeCuts(0),
104 fVtx1(0),
105
106 fHistoDistanceLcToPrimVtx(0),
107 fHistoDistanceV0ToPrimVtx(0),
108 fHistoDistanceV0ToLc(0),
109
110 fHistoDistanceLcToPrimVtxSgn(0),
111 fHistoDistanceV0ToPrimVtxSgn(0),
112 fHistoDistanceV0ToLcSgn(0),
113
114 fHistoVtxLcResidualToPrimVtx(0),
115 fHistoVtxV0ResidualToPrimVtx(0),
116 fHistoVtxV0ResidualToLc(0),
117
118 fHistoMassV0All(0),
119 fHistoDecayLengthV0All(0),
120 fHistoLifeTimeV0All(0),
121
122 fHistoMassV0True(0),
123 fHistoDecayLengthV0True(0),
124 fHistoLifeTimeV0True(0),
125
126 fHistoMassV0TrueFromAOD(0),
127
128 fHistoMassV0TrueK0S(0),
129 fHistoDecayLengthV0TrueK0S(0),
130 fHistoLifeTimeV0TrueK0S(0),
131
132 fHistoMassV0TrueK0SFromAOD(0),
133
134 fHistoMassLcAll(0),
135 fHistoDecayLengthLcAll(0),
136 fHistoLifeTimeLcAll(0),
137
138 fHistoMassLcTrue(0),
139 fHistoDecayLengthLcTrue(0),
140 fHistoLifeTimeLcTrue(0),
141
142 fHistoMassLcTrueFromAOD(0),
143
144 fHistoMassV0fromLcAll(0),
145 fHistoDecayLengthV0fromLcAll(0),
146 fHistoLifeTimeV0fromLcAll(0),
147
148 fHistoMassV0fromLcTrue(0),
149 fHistoDecayLengthV0fromLcTrue(0),
150 fHistoLifeTimeV0fromLcTrue(0),
151
152 fHistoMassLcSgn(0),
153 fHistoMassLcSgnFromAOD(0),
154 fHistoDecayLengthLcSgn(0),
155 fHistoLifeTimeLcSgn(0),
156
157 fHistoMassV0fromLcSgn(0),
158 fHistoDecayLengthV0fromLcSgn(0),
159 fHistoLifeTimeV0fromLcSgn(0),
160
161 fHistoKF(0),
162 fHistoKFV0(0),
163 fHistoKFLc(0),
164
165 fHistoMassKFV0(0),
166 fHistoDecayLengthKFV0(0),
167 fHistoLifeTimeKFV0(0),
168
169 fHistoMassKFLc(0),
170 fHistoDecayLengthKFLc(0),
171 fHistoLifeTimeKFLc(0),
172
173 fHistoArmenterosPodolanskiV0KF(0),
174 fHistoArmenterosPodolanskiV0KFSgn(0),
175 fHistoArmenterosPodolanskiV0AOD(0),
176 fHistoArmenterosPodolanskiV0AODSgn(0),
177
178 fOutputKF(0),
179 fmcLabelLc(-1),
180 ftopoConstraint(kTRUE),
181 fCallKFVertexing(kFALSE),
182 fKeepingOnlyHIJINGBkg(kFALSE),
183 fUtils(0),
ef87f0d6 184 fHistoBackground(0),
185 fCutKFChi2NDF(999999.),
186 fCutKFDeviationFromVtx(999999.),
187 fCutKFDeviationFromVtxV0(0.),
188 fCurrentEvent(-1),
ae48be26 189 fBField(0),
4539d15b 190 fKeepingOnlyPYTHIABkg(kFALSE),
f71760b8 191 fHistoMCLcK0SpGen(0x0),
192 fHistoMCLcK0SpGenAcc(0x0),
294d9511 193 fHistoMCLcK0SpGenLimAcc(0x0),
194 fTriggerMask(0)
0abe5247 195{
196 //
197 // Default ctor
198 //
199}
200//___________________________________________________________________________
201AliAnalysisTaskSELc2V0bachelorTMVA::AliAnalysisTaskSELc2V0bachelorTMVA(const Char_t* name,
202 AliRDHFCutsLctoV0* analCuts, Bool_t useOnTheFly) :
203 AliAnalysisTaskSE(name),
204 fUseMCInfo(kFALSE),
205 fOutput(0),
206 fCEvents(0),
207 fPIDResponse(0),
208 fPIDCombined(0),
209 fIsK0sAnalysis(kFALSE),
210 fCounter(0),
211 fAnalCuts(analCuts),
212 fListCuts(0),
213 fUseOnTheFlyV0(useOnTheFly),
214 fIsEventSelected(kFALSE),
215 fVariablesTreeSgn(0),
216 fVariablesTreeBkg(0),
217 fCandidateVariables(),
218 fIspA(kFALSE),
219 fHistoEvents(0),
220 fHistoLc(0),
221 fHistoLcOnTheFly(0),
222 fFillOnlySgn(kFALSE),
223 fHistoLcBeforeCuts(0),
224 fHistoFiducialAcceptance(0),
225 fHistoCodesSgn(0),
226 fHistoCodesBkg(0),
227 fHistoLcpKpiBeforeCuts(0),
228 fVtx1(0),
229
230 fHistoDistanceLcToPrimVtx(0),
231 fHistoDistanceV0ToPrimVtx(0),
232 fHistoDistanceV0ToLc(0),
233
234 fHistoDistanceLcToPrimVtxSgn(0),
235 fHistoDistanceV0ToPrimVtxSgn(0),
236 fHistoDistanceV0ToLcSgn(0),
237
238 fHistoVtxLcResidualToPrimVtx(0),
239 fHistoVtxV0ResidualToPrimVtx(0),
240 fHistoVtxV0ResidualToLc(0),
241
242 fHistoMassV0All(0),
243 fHistoDecayLengthV0All(0),
244 fHistoLifeTimeV0All(0),
245
246 fHistoMassV0True(0),
247 fHistoDecayLengthV0True(0),
248 fHistoLifeTimeV0True(0),
249
250 fHistoMassV0TrueFromAOD(0),
251
252 fHistoMassV0TrueK0S(0),
253 fHistoDecayLengthV0TrueK0S(0),
254 fHistoLifeTimeV0TrueK0S(0),
255
256 fHistoMassV0TrueK0SFromAOD(0),
257
258 fHistoMassLcAll(0),
259 fHistoDecayLengthLcAll(0),
260 fHistoLifeTimeLcAll(0),
261
262 fHistoMassLcTrue(0),
263 fHistoDecayLengthLcTrue(0),
264 fHistoLifeTimeLcTrue(0),
265
266 fHistoMassLcTrueFromAOD(0),
267
268 fHistoMassV0fromLcAll(0),
269 fHistoDecayLengthV0fromLcAll(0),
270 fHistoLifeTimeV0fromLcAll(0),
271
272 fHistoMassV0fromLcTrue(0),
273 fHistoDecayLengthV0fromLcTrue(0),
274 fHistoLifeTimeV0fromLcTrue(0),
275
276 fHistoMassLcSgn(0),
277 fHistoMassLcSgnFromAOD(0),
278 fHistoDecayLengthLcSgn(0),
279 fHistoLifeTimeLcSgn(0),
280
281 fHistoMassV0fromLcSgn(0),
282 fHistoDecayLengthV0fromLcSgn(0),
283 fHistoLifeTimeV0fromLcSgn(0),
284
285 fHistoKF(0),
286 fHistoKFV0(0),
287 fHistoKFLc(0),
288
289 fHistoMassKFV0(0),
290 fHistoDecayLengthKFV0(0),
291 fHistoLifeTimeKFV0(0),
292
293 fHistoMassKFLc(0),
294 fHistoDecayLengthKFLc(0),
295 fHistoLifeTimeKFLc(0),
296
297 fHistoArmenterosPodolanskiV0KF(0),
298 fHistoArmenterosPodolanskiV0KFSgn(0),
299 fHistoArmenterosPodolanskiV0AOD(0),
300 fHistoArmenterosPodolanskiV0AODSgn(0),
301
302 fOutputKF(0),
303 fmcLabelLc(-1),
304 ftopoConstraint(kTRUE),
305 fCallKFVertexing(kFALSE),
306 fKeepingOnlyHIJINGBkg(kFALSE),
307 fUtils(0),
ef87f0d6 308 fHistoBackground(0),
309 fCutKFChi2NDF(999999.),
310 fCutKFDeviationFromVtx(999999.),
311 fCutKFDeviationFromVtxV0(0.),
312 fCurrentEvent(-1),
ae48be26 313 fBField(0),
4539d15b 314 fKeepingOnlyPYTHIABkg(kFALSE),
f71760b8 315 fHistoMCLcK0SpGen(0x0),
316 fHistoMCLcK0SpGenAcc(0x0),
294d9511 317 fHistoMCLcK0SpGenLimAcc(0x0),
318 fTriggerMask(0)
0abe5247 319{
320 //
321 // Constructor. Initialization of Inputs and Outputs
322 //
323 Info("AliAnalysisTaskSELc2V0bachelorTMVA","Calling Constructor");
324
325 DefineOutput(1, TList::Class()); // Tree signal + Tree Bkg + histoEvents
326 DefineOutput(2, AliNormalizationCounter::Class()); // normalization counter object
327 DefineOutput(3, TList::Class()); // Cuts
328 DefineOutput(4, TTree::Class()); // Tree signal + Tree Bkg + histoEvents
329 DefineOutput(5, TTree::Class()); // Tree signal + Tree Bkg + histoEvents
330 DefineOutput(6, TList::Class()); // Tree signal + Tree Bkg + histoEvents
331
332}
333
334//___________________________________________________________________________
335AliAnalysisTaskSELc2V0bachelorTMVA::~AliAnalysisTaskSELc2V0bachelorTMVA() {
336 //
337 // destructor
338 //
339 Info("~AliAnalysisTaskSELc2V0bachelorTMVA","Calling Destructor");
340
341 if (fOutput) {
342 delete fOutput;
343 fOutput = 0;
344 }
345
346 if (fPIDResponse) {
347 delete fPIDResponse;
348 }
349
350 if (fPIDCombined) {
351 delete fPIDCombined;
352 }
353
354 if (fCounter) {
355 delete fCounter;
356 fCounter = 0;
357 }
358
359 if (fAnalCuts) {
360 delete fAnalCuts;
361 fAnalCuts = 0;
362 }
363
364 if (fListCuts) {
365 delete fListCuts;
366 fListCuts = 0;
367 }
368
369 if(fVariablesTreeSgn){
370 delete fVariablesTreeSgn;
371 fVariablesTreeSgn = 0;
372 }
373
374 if(fVariablesTreeBkg){
375 delete fVariablesTreeBkg;
376 fVariablesTreeBkg = 0;
377 }
378
379 if (fOutputKF) {
380 delete fOutputKF;
381 fOutputKF = 0;
382 }
383
384 if (fUtils) {
385 delete fUtils;
386 fUtils = 0;
387 }
388
389}
390//_________________________________________________
391void AliAnalysisTaskSELc2V0bachelorTMVA::Init() {
392 //
393 // Initialization
394 //
395
396 fIsEventSelected=kFALSE;
397
398 if (fDebug > 1) AliInfo("Init");
399
400 fListCuts = new TList();
401 fListCuts->SetOwner();
402 fListCuts->Add(new AliRDHFCutsLctoV0(*fAnalCuts));
403 PostData(3,fListCuts);
404
ae48be26 405 if (fUseMCInfo && (fKeepingOnlyHIJINGBkg || fKeepingOnlyPYTHIABkg)) fUtils = new AliVertexingHFUtils();
0abe5247 406
407 return;
408}
409
410//________________________________________ terminate ___________________________
411void AliAnalysisTaskSELc2V0bachelorTMVA::Terminate(Option_t*)
412{
413 // The Terminate() function is the last function to be called during
414 // a query. It always runs on the client, it can be used to present
415 // the results graphically or save the results to file.
416
44c36327 417 AliInfo("Terminate");
0abe5247 418 AliAnalysisTaskSE::Terminate();
419
420 fOutput = dynamic_cast<TList*> (GetOutputData(1));
421 if (!fOutput) {
422 AliError("fOutput not available");
423 return;
424 }
4539d15b 425
426
a84144de 427 //AliDebug(2, Form("At MC level, %f Lc --> K0S + p were found", fHistoMCLcK0SpGen->GetEntries()));
428 //AliDebug(2, Form("At MC level, %f Lc --> K0S + p were found in the acceptance", fHistoMCLcK0SpGenAcc->GetEntries()));
429 //AliDebug(2, Form("At Reco level, %lld Lc --> K0S + p were found", fVariablesTreeSgn->GetEntries()));
430 if(fHistoMCLcK0SpGen) {
431 AliInfo(Form("At MC level, %f Lc --> K0S + p were found", fHistoMCLcK0SpGen->GetEntries()));
432 } else {
433 AliInfo("fHistoMCLcK0SpGen not available");
434 }
435 if(fHistoMCLcK0SpGenAcc) {
436 AliInfo(Form("At MC level, %f Lc --> K0S + p were found in the acceptance", fHistoMCLcK0SpGenAcc->GetEntries()));
437 } else {
438 AliInfo("fHistoMCLcK0SpGenAcc not available");
439 }
440 if(fVariablesTreeSgn) {
441 AliInfo(Form("At Reco level, %lld Lc --> K0S + p were found", fVariablesTreeSgn->GetEntries()));
442 } else {
443 AliInfo("fVariablesTreeSgn not available");
444 }
445
0abe5247 446 fOutputKF = dynamic_cast<TList*> (GetOutputData(6));
447 if (!fOutputKF) {
448 AliError("fOutputKF not available");
449 return;
450 }
451
452 return;
453}
454
455//___________________________________________________________________________
456void AliAnalysisTaskSELc2V0bachelorTMVA::UserCreateOutputObjects() {
457 // output
458 AliInfo(Form("CreateOutputObjects of task %s\n", GetName()));
459
460 //slot #1
461 //OpenFile(1);
462 fOutput = new TList();
463 fOutput->SetOwner();
464 fOutput->SetName("listTrees");
465
466 // Output slot 1: list of 2 trees (Sgn + Bkg) of the candidate variables
467 const char* nameoutput = GetOutputSlot(1)->GetContainer()->GetName();
468 fVariablesTreeSgn = new TTree(Form("%s_Sgn", nameoutput), "Candidates variables tree, Signal");
469 fVariablesTreeBkg = new TTree(Form("%s_Bkg", nameoutput), "Candidates variables tree, Background");
ef87f0d6 470 Int_t nVar = 86;
0abe5247 471 fCandidateVariables = new Float_t [nVar];
472 TString * fCandidateVariableNames = new TString[nVar];
473 fCandidateVariableNames[0]="massLc2K0Sp";
474 fCandidateVariableNames[1]="massLc2Lambdapi";
475 fCandidateVariableNames[2]="massK0S";
476 fCandidateVariableNames[3]="massLambda";
477 fCandidateVariableNames[4]="massLambdaBar";
478 fCandidateVariableNames[5]="cosPAK0S";
479 fCandidateVariableNames[6]="dcaV0";
480 fCandidateVariableNames[7]="tImpParBach";
481 fCandidateVariableNames[8]="tImpParV0";
482 fCandidateVariableNames[9]="nSigmaTPCpr";
483 fCandidateVariableNames[10]="nSigmaTPCpi";
484 fCandidateVariableNames[11]="nSigmaTPCka";
485 fCandidateVariableNames[12]="nSigmaTOFpr";
486 fCandidateVariableNames[13]="nSigmaTOFpi";
487 fCandidateVariableNames[14]="nSigmaTOFka";
488 fCandidateVariableNames[15]="bachelorPt";
489 fCandidateVariableNames[16]="V0positivePt";
490 fCandidateVariableNames[17]="V0negativePt";
491 fCandidateVariableNames[18]="dcaV0pos";
492 fCandidateVariableNames[19]="dcaV0neg";
493 fCandidateVariableNames[20]="v0Pt";
494 fCandidateVariableNames[21]="massGamma";
495 fCandidateVariableNames[22]="LcPt";
496 fCandidateVariableNames[23]="combinedProtonProb";
497 fCandidateVariableNames[24]="LcEta";
498 fCandidateVariableNames[25]="V0positiveEta";
499 fCandidateVariableNames[26]="V0negativeEta";
ae48be26 500 fCandidateVariableNames[27]="TPCProtonProb";
501 fCandidateVariableNames[28]="TOFProtonProb";
0abe5247 502 fCandidateVariableNames[29]="bachelorEta";
503 fCandidateVariableNames[30]="LcP";
504 fCandidateVariableNames[31]="bachelorP";
505 fCandidateVariableNames[32]="v0P";
506 fCandidateVariableNames[33]="V0positiveP";
507 fCandidateVariableNames[34]="V0negativeP";
508 fCandidateVariableNames[35]="LcY";
509 fCandidateVariableNames[36]="v0Y";
510 fCandidateVariableNames[37]="bachelorY";
511 fCandidateVariableNames[38]="V0positiveY";
512 fCandidateVariableNames[39]="V0negativeY";
513 fCandidateVariableNames[40]="v0Eta";
514 fCandidateVariableNames[41]="DecayLengthLc";
515 fCandidateVariableNames[42]="DecayLengthK0S";
516 fCandidateVariableNames[43]="CtLc";
517 fCandidateVariableNames[44]="CtK0S";
518 fCandidateVariableNames[45]="bachCode";
519 fCandidateVariableNames[46]="k0SCode";
520
521 fCandidateVariableNames[47]="V0KFmass";
522 fCandidateVariableNames[48]="V0KFdecayLength";
523 fCandidateVariableNames[49]="V0KFlifeTime";
524
525 fCandidateVariableNames[50]="V0KFmassErr";
526 fCandidateVariableNames[51]="V0KFdecayTimeErr";
527 fCandidateVariableNames[52]="V0KFlifeTimeErr";
528
529 fCandidateVariableNames[53]="LcKFmass";
530 fCandidateVariableNames[54]="LcKFdecayLength";
531 fCandidateVariableNames[55]="LcKFlifeTime";
532
533 fCandidateVariableNames[56]="LcKFmassErr";
534 fCandidateVariableNames[57]="LcKFdecayTimeErr";
535 fCandidateVariableNames[58]="LcKFlifeTimeErr";
536
537 fCandidateVariableNames[59]="LcKFDistToPrimVtx";
538 fCandidateVariableNames[60]="V0KFDistToPrimVtx";
539 fCandidateVariableNames[61]="V0KFDistToLc";
540 fCandidateVariableNames[62]="alphaArmKF";
541 fCandidateVariableNames[63]="ptArmKF";
542 fCandidateVariableNames[64]="alphaArm";
543 fCandidateVariableNames[65]="ptArm";
544
ef87f0d6 545 fCandidateVariableNames[66]="ITSrefitV0pos";
546 fCandidateVariableNames[67]="ITSrefitV0neg";
547
548 fCandidateVariableNames[68]="TPCClV0pos";
549 fCandidateVariableNames[69]="TPCClV0neg";
550
551 fCandidateVariableNames[70]="v0Xcoord";
552 fCandidateVariableNames[71]="v0Ycoord";
553 fCandidateVariableNames[72]="v0Zcoord";
554 fCandidateVariableNames[73]="primVtxX";
555 fCandidateVariableNames[74]="primVtxY";
556 fCandidateVariableNames[75]="primVtxZ";
557
558 fCandidateVariableNames[76]="ITSclBach";
559 fCandidateVariableNames[77]="SPDclBach";
560
561 fCandidateVariableNames[78]="ITSclV0pos";
562 fCandidateVariableNames[79]="SPDclV0pos";
563 fCandidateVariableNames[80]="ITSclV0neg";
564 fCandidateVariableNames[81]="SPDclV0neg";
565
566 fCandidateVariableNames[82]="alphaArmLc";
567 fCandidateVariableNames[83]="alphaArmLcCharge";
568 fCandidateVariableNames[84]="ptArmLc";
569
570 fCandidateVariableNames[85]="CosThetaStar";
571
572
0abe5247 573 for(Int_t ivar=0; ivar<nVar; ivar++){
574 fVariablesTreeSgn->Branch(fCandidateVariableNames[ivar].Data(),&fCandidateVariables[ivar],Form("%s/f",fCandidateVariableNames[ivar].Data()));
575 fVariablesTreeBkg->Branch(fCandidateVariableNames[ivar].Data(),&fCandidateVariables[ivar],Form("%s/f",fCandidateVariableNames[ivar].Data()));
576 }
577
578 fHistoEvents = new TH1F("fHistoEvents", "fHistoEvents", 2, -0.5, 1.5);
579 TString labelEv[2] = {"NotSelected", "Selected"};
580 for (Int_t ibin = 1; ibin <= fHistoEvents->GetNbinsX(); ibin++){
581 fHistoEvents->GetXaxis()->SetBinLabel(ibin, labelEv[ibin-1].Data());
582 }
583
584 fHistoLc = new TH1F("fHistoLc", "fHistoLc", 2, -0.5, 1.5);
585
586 fHistoLcOnTheFly = new TH1F("fHistoLcOnTheFly", "fHistoLcOnTheFly", 4, -0.5, 3.5);
587 TString labelOnTheFly[4] = {"OnTheFly-Bkg", "OfflineBkg", "OnTheFly-Sgn", "OfflineSgn"};
588 for (Int_t ibin = 1; ibin <= fHistoLcOnTheFly->GetNbinsX(); ibin++){
589 fHistoLcOnTheFly->GetXaxis()->SetBinLabel(ibin, labelOnTheFly[ibin-1].Data());
590 }
591
592 fHistoLcBeforeCuts = new TH1F("fHistoLcBeforeCuts", "fHistoLcBeforeCuts", 2, -0.5, 1.5);
593 TString labelBeforeCuts[2] = {"Bkg", "Sgn"};
594 for (Int_t ibin = 1; ibin <= fHistoLcBeforeCuts->GetNbinsX(); ibin++){
595 fHistoLcBeforeCuts->GetXaxis()->SetBinLabel(ibin, labelBeforeCuts[ibin-1].Data());
596 fHistoLc->GetXaxis()->SetBinLabel(ibin, labelBeforeCuts[ibin-1].Data());
597 }
598
599 fHistoFiducialAcceptance = new TH1F("fHistoFiducialAcceptance", "fHistoFiducialAcceptance", 4, -0.5, 3.5);
600 TString labelAcc[4] = {"NotOk-Bkg", "Ok-Bkg", "NotOk-Sgn", "Ok-Sgn"};
601 for (Int_t ibin = 1; ibin <= fHistoFiducialAcceptance->GetNbinsX(); ibin++){
602 fHistoFiducialAcceptance->GetXaxis()->SetBinLabel(ibin, labelAcc[ibin-1].Data());
603 }
604
605 fHistoCodesSgn = new TH2F("fHistoCodesSgn", "fHistoCodes for Signal; bachelor; K0S", 7, -1.5, 5.5, 9, -1.5, 7.5);
606 fHistoCodesBkg = new TH2F("fHistoCodesBkg", "fHistoCodes for Background; bachelor; K0S", 7, -1.5, 5.5, 9, -1.5, 7.5);
607
608 TString labelx[7] = { "kBachInvalid", "kBachFake", "kBachNoProton", "kBachPrimary", "kBachNoLambdaMother",
609 "kBachDifferentLambdaMother", "kBachCorrectLambdaMother"};
610 TString labely[9] = { "kK0SInvalid", "kK0SFake", "kK0SNoK0S", "kK0SWithoutMother", "kK0SNotFromK0",
611 "kK0Primary", "kK0NoLambdaMother", "kK0DifferentLambdaMother", "kK0CorrectLambdaMother"};
612
613 for (Int_t ibin = 1; ibin <= fHistoCodesSgn->GetNbinsX(); ibin++){
614 fHistoCodesSgn->GetXaxis()->SetBinLabel(ibin, labelx[ibin-1].Data());
615 fHistoCodesBkg->GetXaxis()->SetBinLabel(ibin, labelx[ibin-1].Data());
616 }
617 for (Int_t ibin = 1; ibin <= fHistoCodesSgn->GetNbinsY(); ibin++){
618 fHistoCodesSgn->GetYaxis()->SetBinLabel(ibin, labely[ibin-1].Data());
619 fHistoCodesBkg->GetYaxis()->SetBinLabel(ibin, labely[ibin-1].Data());
620 }
621
622 fHistoLcpKpiBeforeCuts = new TH1F("fHistoLcpKpiBeforeCuts", "fHistoLcpKpiBeforeCuts", 2, -0.5, 1.5);
623 for (Int_t ibin = 1; ibin <= fHistoLcpKpiBeforeCuts->GetNbinsX(); ibin++){
624 fHistoLcpKpiBeforeCuts->GetXaxis()->SetBinLabel(ibin, labelBeforeCuts[ibin-1].Data());
625 }
626
ae48be26 627 fHistoBackground = new TH1F("fHistoBackground", "fHistoBackground", 4, -0.5, 3.5);
628 TString labelBkg[4] = {"Injected", "Non-injected", "Non-PYTHIA", "PYTHIA"};
0abe5247 629 for (Int_t ibin = 1; ibin <= fHistoBackground->GetNbinsX(); ibin++){
630 fHistoBackground->GetXaxis()->SetBinLabel(ibin, labelBkg[ibin-1].Data());
631 }
632
633 //fOutput->Add(fVariablesTreeSgn);
634 //fOutput->Add(fVariablesTreeBkg);
635
4539d15b 636 const Float_t ptbins[15] = {0., 1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 12., 17., 25., 35.};
637
f71760b8 638 fHistoMCLcK0SpGen = new TH1F("fHistoMCLcK0SpGen", "fHistoMCLcK0SpGen", 14, ptbins);
639 fHistoMCLcK0SpGenAcc = new TH1F("fHistoMCLcK0SpGenAcc", "fHistoMCLcK0SpGenAcc", 14, ptbins);
640 fHistoMCLcK0SpGenLimAcc = new TH1F("fHistoMCLcK0SpGenLimAcc", "fHistoMCLcK0SpGenLimAcc", 14, ptbins);
4539d15b 641
0abe5247 642 fOutput->Add(fHistoEvents);
643 fOutput->Add(fHistoLc);
644 fOutput->Add(fHistoLcOnTheFly);
645 fOutput->Add(fHistoLcBeforeCuts);
646 fOutput->Add(fHistoFiducialAcceptance);
647 fOutput->Add(fHistoCodesSgn);
648 fOutput->Add(fHistoCodesBkg);
649 fOutput->Add(fHistoLcpKpiBeforeCuts);
650 fOutput->Add(fHistoBackground);
f71760b8 651 fOutput->Add(fHistoMCLcK0SpGen);
652 fOutput->Add(fHistoMCLcK0SpGenAcc);
653 fOutput->Add(fHistoMCLcK0SpGenLimAcc);
0abe5247 654
655 PostData(1, fOutput);
656 PostData(4, fVariablesTreeSgn);
657 PostData(5, fVariablesTreeBkg);
658 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
659 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
660 fPIDResponse = inputHandler->GetPIDResponse();
661
662 if (fAnalCuts->GetIsUsePID()){
663 /*
664 fAnalCuts->GetPidHF()->SetPidResponse(fPIDResponse);
665 fAnalCuts->GetPidV0pos()->SetPidResponse(fPIDResponse);
666 fAnalCuts->GetPidV0neg()->SetPidResponse(fPIDResponse);
667 fAnalCuts->GetPidHF()->SetOldPid(kFALSE);
668 fAnalCuts->GetPidV0pos()->SetOldPid(kFALSE);
669 fAnalCuts->GetPidV0neg()->SetOldPid(kFALSE);
670 */
671 fAnalCuts->SetUsePID(kFALSE); // I don't want to use the PID through the cut object, but I will use the PID response directly!!!
672 }
673
674 // Setting properties of PID
675 fPIDCombined=new AliPIDCombined;
676 fPIDCombined->SetDefaultTPCPriors();
677 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC+AliPIDResponse::kDetTOF);
678 //fPIDCombined->SetPriorDistribution((AliPID::EParticleType)ispec,fPriors[ispec]);
679
0abe5247 680 fCounter = new AliNormalizationCounter("NormalizationCounter");
681 fCounter->Init();
682 PostData(2, fCounter);
683
684 // Histograms from KF
685
ef87f0d6 686 if (fCallKFVertexing){
687 fHistoDistanceLcToPrimVtx = new TH1D("fHistoDistanceLcToPrimVtx", "Lc distance to Prim Vertex from KF; distance [cm]", 1000, 0., 1);
688 fHistoDistanceV0ToPrimVtx = new TH1D("fHistoDistanceV0ToPrimVtx", "V0 distance to Prim Vertex from KF; distance [cm]", 1000, 0., 100.);
689 fHistoDistanceV0ToLc = new TH1D("fHistoDistanceV0ToLc", "V0 distance to Lc Vertex from KF; distance [cm]", 1000, 0., 100.);
690
691 fHistoDistanceLcToPrimVtxSgn = new TH1D("fHistoDistanceLcToPrimVtxSgn", "Lc Sgn distance to Prim Vertex from KF; distance [cm]", 1000, 0., 1);
692 fHistoDistanceV0ToPrimVtxSgn = new TH1D("fHistoDistanceV0ToPrimVtxSgn", "V0 Sgn distance to Prim Vertex from KF; distance [cm]", 1000, 0., 100.);
693 fHistoDistanceV0ToLcSgn = new TH1D("fHistoDistanceV0ToLcSgn", "V0 Sgn distance to Lc Vertex from KF; distance [cm]", 1000, 0., 100.);
694
695 fHistoVtxLcResidualToPrimVtx = new TH1D("fHistoVtxLcResidualToPrimVtx", "Residual between MC and KF (MC - KF): Lc to Prim Vtx; distance [cm]", 1000, -5., 5.);
696 fHistoVtxV0ResidualToPrimVtx = new TH1D("fHistoVtxV0ResidualToPrimVtx", "Residual between MC and KF (MC - KF): V0 to Prim Vtx; distance [cm]", 1000, -5., 5.);
697 fHistoVtxV0ResidualToLc = new TH1D("fHistoVtxV0ResidualToLc", "Residual between MC and KF: V0 to Lc (MC - KF); distance [cm]", 1000, -5., 5.);
698
699 fHistoMassV0All = new TH1D("fHistoMassV0All", "V0 Mass; mass", 500, 0.4, 0.6);
700 fHistoDecayLengthV0All = new TH1D("fHistoDecayLengthV0All", "V0 Decay Length; decayLength", 500, -10, 10.0);
701 fHistoLifeTimeV0All = new TH1D("fHistoLifeTimeV0All", "V0 Life Time; lifeTime", 500, -10.0, 10.0);
702
703 fHistoMassV0True = new TH1D("fHistoMassV0True", "True V0 Mass; mass", 500, 0.4, 0.6);
704 fHistoDecayLengthV0True = new TH1D("fHistoDecayLengthV0True", "True V0 Decay Length; decayLength", 500, -10, 10.0);
705 fHistoLifeTimeV0True = new TH1D("fHistoLifeTimeV0True", "True V0 Life Time; lifeTime", 500, -10.0, 10.0);
706
707 fHistoMassV0TrueFromAOD = new TH1D("fHistoMassV0TrueFormAOD", "True V0 Mass (AOD); mass", 500, 0.4, 0.6);
708
709 fHistoMassV0TrueK0S = new TH1D("fHistoMassV0TrueK0S", "True V0-K0S Mass; mass", 500, 0.4, 0.6);
710 fHistoDecayLengthV0TrueK0S = new TH1D("fHistoDecayLengthV0TrueK0S", "True V0-K0S Decay Length; decayLength", 500, -10, 10.0);
711 fHistoLifeTimeV0TrueK0S = new TH1D("fHistoLifeTimeV0TrueK0S", "True V0-K0S Life Time; lifeTime", 500, -10.0, 10.0);
712
713 fHistoMassV0TrueK0SFromAOD = new TH1D("fHistoMassV0TrueK0SFormAOD", "True V0-K0S Mass (AOD); mass", 500, 0.4, 0.6);
714
715 fHistoMassLcAll = new TH1D("fHistoMassLcAll", "Lc Mass; mass", 500, 2.0, 3.0);
716 fHistoDecayLengthLcAll = new TH1D("fHistoDecayLengthLcAll", "Lc Decay Length; decayLenght", 100000, -0.1, 0.1);
717 fHistoLifeTimeLcAll = new TH1D("fHistoLifeTimeLcAll", "Lc Life Time; lifeTime", 100000, -0.1, 0.1);
718
719 fHistoMassLcTrue = new TH1D("fHistoMassLcTrue", "True Lc Mass; mass", 500, 2.0, 3.0);
720 fHistoDecayLengthLcTrue = new TH1D("fHistoDecayLengthLcTrue", "True Lc Decay Length; decayLength", 100000, -0.1, 0.1);
721 fHistoLifeTimeLcTrue = new TH1D("fHistoLifeTimeLcTrue", "True Lc Life Time; lifeTime", 100000, -0.1, 0.1);
722
723 fHistoMassLcTrueFromAOD = new TH1D("fHistoMassLcTrueFromAOD", "True Lc Mass (AOD); mass", 500, 2.0, 3.0);
724
725 fHistoMassV0fromLcAll = new TH1D("fHistoMassV0fromLcAll", "V0 mass from Lc built in KF; mass", 500, 0.4, 0.6);
726 fHistoDecayLengthV0fromLcAll = new TH1D("fHistoDecayLengthV0fromLcAll", "V0 Decay Length from Lc built in KF; decayLength", 500, 0, 10.0);
727 fHistoLifeTimeV0fromLcAll = new TH1D("fHistoLifeTimeV0fromLcAll", "V0 Life Time from Lc built in KF; lifeTime", 500, 0.0, 3.0);
728
729 fHistoMassV0fromLcTrue = new TH1D("fHistoMassV0fromLcTrue", "V0 mass from true Lc built in KF; mass", 500, 0.4, 0.6);
730 fHistoDecayLengthV0fromLcTrue= new TH1D("fHistoDecayLengthV0fromLcTrue", "V0 Decay Length from true Lc built in KF; decayLength", 500, 0, 10.0);
731 fHistoLifeTimeV0fromLcTrue = new TH1D("fHistoLifeTimeV0fromLcTrue", "V0 Life Time from true Lc built in KF; lifeTime", 500, 0.0, 3.0);
732
733 fHistoMassLcSgn = new TH1D("fHistoMassLcSgn", "True Lc Signal Mass; mass", 500, 2.0, 3.0);
734 fHistoMassLcSgnFromAOD = new TH1D("fHistoMassLcSgnFromAOD", "True Lc Signal Mass (AOD); mass", 500, 2.0, 3.0);
735 fHistoDecayLengthLcSgn = new TH1D("fHistoDecayLengthLcSgn", "True Lc Signal Decay Length; decayLength", 100000, -0.1, 0.1);
736 fHistoLifeTimeLcSgn = new TH1D("fHistoLifeTimeLcSgn", "True Lc Signal Life Time; lifeTime", 100000, -0.1, 0.1);
737
738 fHistoMassV0fromLcSgn = new TH1D("fHistoMassV0fromLcSgn", "V0 from True Lc Signal Mass; mass", 500, 0.4, 0.6);
739 fHistoDecayLengthV0fromLcSgn = new TH1D("fHistoDecayLengthV0fromLcSgn", "V0 True Lc Signal Decay Length; decayLength", 500, 0, 10.0);
740 fHistoLifeTimeV0fromLcSgn = new TH1D("fHistoLifeTimeV0fromLcSgn", "V0 True Lc Signal Life Time; lifeTime", 500, 0.0, 3.0);
741
742 fHistoKF = new TH2D("fHistoKF", "Summary from KF; V0 KF; Lc KF", 16, -0.5, 15.5, 16, -0.5, 15.5);
743 fHistoKFV0 = new TH1D("fHistoKFV0", "Summary from KF; V0 KF", 16, -0.5, 15.5);
744 fHistoKFLc = new TH1D("fHistoKFLc", "Summary from KF; V0 KF", 16, -0.5, 15.5);
745 TString axisLabel[16] = {"AllOk", "M_NotOk", "Sm_NotOk", "Dl_NotOk",
746 "Lt_NotOk", "M_Sm_NotOk", "M_Dl_NotOk", "M_Lt_NotOk",
747 "Dl_Sm_NotOk", "Dl_Lt_NotOk", "Sm_Lt_NotOk", "M_Sm_Dl_NotOk",
748 "M_Sm_Lt_NotOk", "Sm_Dl_Lt_NotOk", "M_Dl_Lt_NotOk", "All_NotOk"};
749
750 for (Int_t ibin = 1; ibin <=16; ibin++){
751 fHistoKF->GetXaxis()->SetBinLabel(ibin, axisLabel[ibin-1].Data());
752 fHistoKF->GetYaxis()->SetBinLabel(ibin, axisLabel[ibin-1].Data());
753 fHistoKFV0->GetXaxis()->SetBinLabel(ibin, axisLabel[ibin-1].Data());
754 fHistoKFLc->GetXaxis()->SetBinLabel(ibin, axisLabel[ibin-1].Data());
755 }
756
757 fHistoMassKFV0 = new TH2D("fHistoMassKFV0", "mass vs sigmaMass for V0; mass; sigmaMass", 500, 0.4, 0.6, 500, 0., 10);
758 fHistoDecayLengthKFV0 = new TH2D("fHistoDecayLengthKFV0", "decayLength vs sigmaDecayLength for V0; decayLength; sigmaDecayLength", 500, -10, 10, 500, 0., 10);
759 fHistoLifeTimeKFV0 = new TH2D("fHistoLifeTimeKFV0", "lifeTime vs sigmalifeTime for V0; lifeTime; sigmaLifeTime", 500, -10, 10, 500, 0., 10);
760
761 fHistoMassKFLc = new TH2D("fHistoMassKFLc", "mass vs sigmaMass for Lc; mass; sigmaMass", 500, 0.4, 0.6, 500, 0., 10);
762 fHistoDecayLengthKFLc = new TH2D("fHistoDecayLengthKFLc", "decayLength vs sigmaDecayLength for Lc; decayLength; sigmaDecayLength", 500, -10, 10, 500, 0., 10);
763 fHistoLifeTimeKFLc = new TH2D("fHistoLifeTimeKFLc", "lifeTime vs sigmalifeTime for Lc; lifeTime; sigmaLifeTime", 500, -10, 10, 500, 0., 10);
764
765 fHistoArmenterosPodolanskiV0KF = new TH2D("fHistoArmenterosPodolanskiV0KF", "V0 ArmenterosPodolanski from KF; #alpha; Qt", 1000, -1, 1, 1000, 0, 1);
766 fHistoArmenterosPodolanskiV0KFSgn = new TH2D("fHistoArmenterosPodolanskiV0KFSgn", "V0 (signal) ArmenterosPodolanski from KF; #alpha; Qt", 1000, -1, 1, 1000, 0, 1);
767 fHistoArmenterosPodolanskiV0AOD = new TH2D("fHistoArmenterosPodolanskiV0AOD", "V0 ArmenterosPodolanski from AOD; #alpha; Qt", 1000, -1, 1, 1000, 0, 1);
768 fHistoArmenterosPodolanskiV0AODSgn = new TH2D("fHistoArmenterosPodolanskiV0AODSgn", "V0 (signal) ArmenterosPodolanski from AOD; #alpha; Qt", 1000, -1, 1, 1000, 0, 1);
0abe5247 769 }
770
0abe5247 771 fOutputKF = new TList();
772 fOutputKF->SetOwner();
773 fOutputKF->SetName("listHistoKF");
774
ef87f0d6 775 if (fCallKFVertexing){
776 fOutputKF->Add(fHistoDistanceLcToPrimVtx);
777 fOutputKF->Add(fHistoDistanceV0ToPrimVtx);
778 fOutputKF->Add(fHistoDistanceV0ToLc);
779
780 fOutputKF->Add(fHistoDistanceLcToPrimVtxSgn);
781 fOutputKF->Add(fHistoDistanceV0ToPrimVtxSgn);
782 fOutputKF->Add(fHistoDistanceV0ToLcSgn);
783
784 fOutputKF->Add(fHistoVtxLcResidualToPrimVtx);
785 fOutputKF->Add(fHistoVtxV0ResidualToPrimVtx);
786 fOutputKF->Add(fHistoVtxV0ResidualToLc);
787
788 fOutputKF->Add(fHistoMassV0All);
789 fOutputKF->Add(fHistoDecayLengthV0All);
790 fOutputKF->Add(fHistoLifeTimeV0All);
791
792 fOutputKF->Add(fHistoMassV0True);
793 fOutputKF->Add(fHistoDecayLengthV0True);
794 fOutputKF->Add(fHistoLifeTimeV0True);
795
796 fOutputKF->Add(fHistoMassV0TrueFromAOD);
797
798 fOutputKF->Add(fHistoMassV0TrueK0S);
799 fOutputKF->Add(fHistoDecayLengthV0TrueK0S);
800 fOutputKF->Add(fHistoLifeTimeV0TrueK0S);
801
802 fOutputKF->Add(fHistoMassV0TrueK0SFromAOD);
803
804 fOutputKF->Add(fHistoMassLcAll);
805 fOutputKF->Add(fHistoDecayLengthLcAll);
806 fOutputKF->Add(fHistoLifeTimeLcAll);
807
808 fOutputKF->Add(fHistoMassLcTrue);
809 fOutputKF->Add(fHistoDecayLengthLcTrue);
810 fOutputKF->Add(fHistoLifeTimeLcTrue);
811
812 fOutputKF->Add(fHistoMassLcTrueFromAOD);
813
814 fOutputKF->Add(fHistoMassV0fromLcAll);
815 fOutputKF->Add(fHistoDecayLengthV0fromLcAll);
816 fOutputKF->Add(fHistoLifeTimeV0fromLcAll);
817
818 fOutputKF->Add(fHistoMassV0fromLcTrue);
819 fOutputKF->Add(fHistoDecayLengthV0fromLcTrue);
820 fOutputKF->Add(fHistoLifeTimeV0fromLcTrue);
821
822 fOutputKF->Add(fHistoMassLcSgn);
823 fOutputKF->Add(fHistoMassLcSgnFromAOD);
824 fOutputKF->Add(fHistoDecayLengthLcSgn);
825 fOutputKF->Add(fHistoLifeTimeLcSgn);
826
827 fOutputKF->Add(fHistoMassV0fromLcSgn);
828 fOutputKF->Add(fHistoDecayLengthV0fromLcSgn);
829 fOutputKF->Add(fHistoLifeTimeV0fromLcSgn);
830
831 fOutputKF->Add(fHistoKF);
832 fOutputKF->Add(fHistoKFV0);
833 fOutputKF->Add(fHistoKFLc);
834
835 fOutputKF->Add(fHistoMassKFV0);
836 fOutputKF->Add(fHistoDecayLengthKFV0);
837 fOutputKF->Add(fHistoLifeTimeKFV0);
838
839 fOutputKF->Add(fHistoMassKFLc);
840 fOutputKF->Add(fHistoDecayLengthKFLc);
841 fOutputKF->Add(fHistoLifeTimeKFLc);
842
843 fOutputKF->Add(fHistoArmenterosPodolanskiV0KF);
844 fOutputKF->Add(fHistoArmenterosPodolanskiV0KFSgn);
845 fOutputKF->Add(fHistoArmenterosPodolanskiV0AOD);
846 fOutputKF->Add(fHistoArmenterosPodolanskiV0AODSgn);
847 }
0abe5247 848
849 PostData(6, fOutputKF);
850
0abe5247 851 return;
852}
853
854//_________________________________________________
855void AliAnalysisTaskSELc2V0bachelorTMVA::UserExec(Option_t *)
856{
857 // user exec
858 if (!fInputEvent) {
859 AliError("NO EVENT FOUND!");
860 return;
861 }
862
ef87f0d6 863 fCurrentEvent++;
4539d15b 864 AliDebug(2, Form("Processing event = %d", fCurrentEvent));
0abe5247 865 AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
866 TClonesArray *arrayLctopKos=0;
867
868 TClonesArray *array3Prong = 0;
869
870 if (!aodEvent && AODEvent() && IsStandardAOD()) {
871 // In case there is an AOD handler writing a standard AOD, use the AOD
872 // event in memory rather than the input (ESD) event.
873 aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());
874 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
875 // have to taken from the AOD event hold by the AliAODExtension
876 AliAODHandler* aodHandler = (AliAODHandler*)
877 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
878
879 if (aodHandler->GetExtensions()) {
880 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
881 AliAODEvent *aodFromExt = ext->GetAOD();
882 arrayLctopKos=(TClonesArray*)aodFromExt->GetList()->FindObject("CascadesHF");
883
884 array3Prong=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
885 }
886 } else {
887 arrayLctopKos=(TClonesArray*)aodEvent->GetList()->FindObject("CascadesHF");
888
889 array3Prong=(TClonesArray*)aodEvent->GetList()->FindObject("Charm3Prong");
890 }
891
892 if ( !fUseMCInfo && fIspA) {
893 fAnalCuts->SetTriggerClass("");
294d9511 894 fAnalCuts->SetTriggerMask(fTriggerMask);
0abe5247 895 }
bd95bc1f 896
897 Int_t runnumber = aodEvent->GetRunNumber();
898 if (aodEvent->GetTriggerMask() == 0 && (runnumber >= 195344 && runnumber <= 195677)){
899 AliDebug(3,"Event rejected because of null trigger mask");
900 return;
901 }
902
0abe5247 903 fCounter->StoreEvent(aodEvent,fAnalCuts,fUseMCInfo);
bd95bc1f 904
0abe5247 905 // mc analysis
906 TClonesArray *mcArray = 0;
907 AliAODMCHeader *mcHeader=0;
908
909 if (fUseMCInfo) {
f71760b8 910 // MC array need for matching
0abe5247 911 mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));
912 if (!mcArray) {
913 AliError("Could not find Monte-Carlo in AOD");
914 return;
915 }
916 // load MC header
917 mcHeader = (AliAODMCHeader*)aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName());
918 if (!mcHeader) {
919 AliError("AliAnalysisTaskSELc2V0bachelorTMVA::UserExec: MC header branch not found!\n");
920 return;
921 }
47a1c02d 922
923 Double_t zMCVertex = mcHeader->GetVtxZ();
924 if (TMath::Abs(zMCVertex) > fAnalCuts->GetMaxVtxZ()){
925 AliDebug(3,Form("z coordinate of MC vertex = %f, it was required to be within [-%f, +%f], skipping event", zMCVertex, fAnalCuts->GetMaxVtxZ(), fAnalCuts->GetMaxVtxZ()));
926 AliInfo(Form("z coordinate of MC vertex = %f, it was required to be within [-%f, +%f], skipping event", zMCVertex, fAnalCuts->GetMaxVtxZ(), fAnalCuts->GetMaxVtxZ()));
927 return;
928 }
929
4539d15b 930 //Printf("Filling MC histo");
931 FillMCHisto(mcArray);
932 }
47a1c02d 933
4539d15b 934 // AOD primary vertex
935 fVtx1 = (AliAODVertex*)aodEvent->GetPrimaryVertex();
936 if (!fVtx1) return;
937 if (fVtx1->GetNContributors()<1) return;
938
939 fIsEventSelected = fAnalCuts->IsEventSelected(aodEvent);
940
941 if ( !fIsEventSelected ) {
942 fHistoEvents->Fill(0);
943 return; // don't take into account not selected events
0abe5247 944 }
4539d15b 945 fHistoEvents->Fill(1);
946
947 // Setting magnetic field for KF vertexing
948 fBField = aodEvent->GetMagneticField();
949 AliKFParticle::SetField(fBField);
0abe5247 950
951 Int_t nSelectedAnal = 0;
952 if (fIsK0sAnalysis) {
953 MakeAnalysisForLc2prK0S(arrayLctopKos, mcArray,
954 nSelectedAnal, fAnalCuts,
955 array3Prong, mcHeader);
956 }
957 fCounter->StoreCandidates(aodEvent,nSelectedAnal,kFALSE);
958
959 PostData(1, fOutput);
960 PostData(2, fCounter);
961 PostData(4, fVariablesTreeSgn);
962 PostData(5, fVariablesTreeBkg);
963 PostData(6, fOutputKF);
964
4539d15b 965}
966//-------------------------------------------------------------------------------
967void AliAnalysisTaskSELc2V0bachelorTMVA::FillMCHisto(TClonesArray *mcArray){
968
969 // method to fill MC histo: how many Lc --> K0S + p are there at MC level
970 for (Int_t iPart=0; iPart<mcArray->GetEntriesFast(); iPart++) {
971 AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(mcArray->At(iPart));
972 if (!mcPart){
973 AliError("Failed casting particle from MC array!, Skipping particle");
4539d15b 974 continue;
975 }
976 Int_t pdg = mcPart->GetPdgCode();
977 if (TMath::Abs(pdg) != 4122){
978 AliDebug(2, Form("MC particle %d is not a Lc: its pdg code is %d", iPart, pdg));
979 continue;
980 }
f3c24752 981 AliDebug(2, Form("Step 0 ok: MC particle %d is a Lc: its pdg code is %d", iPart, pdg));
4539d15b 982 Int_t labeldaugh0 = mcPart->GetDaughter(0);
983 Int_t labeldaugh1 = mcPart->GetDaughter(1);
984 if (labeldaugh0 <= 0 || labeldaugh1 <= 0){
985 AliDebug(2, Form("The MC particle doesn't have correct daughters, skipping!!"));
986 continue;
987 }
988 else if (labeldaugh1 - labeldaugh0 == 1){
f3c24752 989 AliDebug(2, Form("Step 1 ok: The MC particle has correct daughters!!"));
4539d15b 990 AliAODMCParticle* daugh0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labeldaugh0));
991 AliAODMCParticle* daugh1 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labeldaugh1));
61284b3b 992 if(!daugh0 || !daugh1){
993 AliDebug(2,"Particle daughters not properly retrieved!");
994 return;
995 }
4539d15b 996 Int_t pdgCodeDaugh0 = TMath::Abs(daugh0->GetPdgCode());
997 Int_t pdgCodeDaugh1 = TMath::Abs(daugh1->GetPdgCode());
f71760b8 998 AliAODMCParticle* bachelorMC = daugh0;
4539d15b 999 AliAODMCParticle* v0MC = daugh1;
1000 AliDebug(2, Form("pdgCodeDaugh0 = %d, pdgCodeDaugh1 = %d", pdgCodeDaugh0, pdgCodeDaugh1));
1001 if ((pdgCodeDaugh0 == 311 && pdgCodeDaugh1 == 2212) || (pdgCodeDaugh0 == 2212 && pdgCodeDaugh1 == 311)){
1002 // we are in the case of Lc --> K0 + p; now we have to check if the K0 decays in K0S, and if this goes in pi+pi-
1003 /// first, we set the bachelor and the v0: above we assumed first proton and second V0, but we could have to change it:
1004 if (pdgCodeDaugh0 == 311 && pdgCodeDaugh1 == 2212) {
f71760b8 1005 bachelorMC = daugh1;
4539d15b 1006 v0MC = daugh0;
1007 }
1008 AliDebug(2, Form("Number of Daughters of v0 = %d", v0MC->GetNDaughters()));
1009 if (v0MC->GetNDaughters() != 1) {
1010 AliDebug(2, "The K0 does not decay in 1 body only! Impossible... Continuing...");
1011 continue;
1012 }
1013 else { // So far: Lc --> K0 + p, K0 with 1 daughter
f3c24752 1014 AliDebug(2, "Step 2 ok: The K0 does decay in 1 body only! ");
4539d15b 1015 Int_t labelK0daugh = v0MC->GetDaughter(0);
1016 AliAODMCParticle* partK0S = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelK0daugh));
1017 if(!partK0S){
1018 AliError("Error while casting particle! returning a NULL array");
4539d15b 1019 continue;
1020 }
1021 else { // So far: Lc --> K0 + p, K0 with 1 daughter that we can access
1022 if (partK0S->GetNDaughters() != 2 || TMath::Abs(partK0S->GetPdgCode() != 310)){
1023 AliDebug(2, "The K0 daughter is not a K0S or does not decay in 2 bodies");
1024 continue;
1025 }
1026 else { // So far: Lc --> K0 + p, K0 --> K0S, K0S in 2 bodies
f3c24752 1027 AliDebug(2, "Step 3 ok: The K0 daughter is a K0S and does decay in 2 bodies");
4539d15b 1028 Int_t labelK0Sdaugh0 = partK0S->GetDaughter(0);
1029 Int_t labelK0Sdaugh1 = partK0S->GetDaughter(1);
1030 AliAODMCParticle* daughK0S0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelK0Sdaugh0));
1031 AliAODMCParticle* daughK0S1 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelK0Sdaugh1));
1032 if (!daughK0S0 || ! daughK0S1){
1033 AliDebug(2, "Could not access K0S daughters, continuing...");
1034 continue;
1035 }
1036 else { // So far: Lc --> K0 + p, K0 --> K0S, K0S in 2 bodies that we can access
f3c24752 1037 AliDebug(2, "Step 4 ok: Could access K0S daughters, continuing...");
4539d15b 1038 Int_t pdgK0Sdaugh0 = daughK0S0->GetPdgCode();
1039 Int_t pdgK0Sdaugh1 = daughK0S1->GetPdgCode();
1040 if (TMath::Abs(pdgK0Sdaugh0) != 211 || TMath::Abs(pdgK0Sdaugh1) != 211){
1041 AliDebug(2, "The K0S does not decay in pi+pi-, continuing");
f3c24752 1042 //AliInfo("The K0S does not decay in pi+pi-, continuing");
4539d15b 1043 }
1044 else { // Full chain: Lc --> K0 + p, K0 --> K0S, K0S --> pi+pi-
1045 if (fAnalCuts->IsInFiducialAcceptance(mcPart->Pt(), mcPart->Y())) {
1046 AliDebug(2, Form("----> Filling histo with pt = %f", mcPart->Pt()));
f71760b8 1047 if(TMath::Abs(mcPart->Y()) < 0.5) fHistoMCLcK0SpGenLimAcc->Fill(mcPart->Pt());
44c36327 1048 //AliInfo(Form("\nparticle = %d, Filling MC Gen histo\n", iPart));
f71760b8 1049 fHistoMCLcK0SpGen->Fill(mcPart->Pt());
1050 if(!(TMath::Abs(bachelorMC->Eta()) > 0.9 || bachelorMC->Pt() < 0.1 ||
1051 TMath::Abs(daughK0S0->Eta()) > 0.9 || daughK0S0->Pt() < 0.1 ||
1052 TMath::Abs(daughK0S1->Eta()) > 0.9 || daughK0S1->Pt() < 0.1)) {
1053 fHistoMCLcK0SpGenAcc->Fill(mcPart->Pt());
1054 }
4539d15b 1055 }
1056 else {
1057 AliDebug(2, "not in fiducial acceptance! Skipping");
1058 continue;
1059 }
1060 }
1061 }
1062 }
1063 }
1064 }
1065 }
1066 }
1067 } // closing loop over mcArray
1068
1069 return;
1070
0abe5247 1071}
1072
1073//-------------------------------------------------------------------------------
1074void AliAnalysisTaskSELc2V0bachelorTMVA::MakeAnalysisForLc2prK0S(TClonesArray *arrayLctopKos,
1075 TClonesArray *mcArray,
1076 Int_t &nSelectedAnal,
1077 AliRDHFCutsLctoV0 *cutsAnal, TClonesArray *array3Prong,
1078 AliAODMCHeader* aodheader){
1079 //Lc prong needed to MatchToMC method
1080
1081 Int_t pdgCand = 4122;
1082 Int_t pdgDgLctoV0bachelor[2]={2212, 310};
1083 Int_t pdgDgV0toDaughters[2]={211, 211};
1084
1085 Int_t pdgDgLctopKpi[3]={2212, 321, 211};
1086
1087 // loop to search for candidates Lc->p+K+pi
1088 Int_t n3Prong = array3Prong->GetEntriesFast();
1089 Int_t nCascades= arrayLctopKos->GetEntriesFast();
1090
1091 //AliInfo(Form("\n\n\n\n3 prong candidates = %d, ncascades = %d \n\n\n\n\n", n3Prong, nCascades));
1092 for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
1093 AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
1094 //Filling a control histogram with no cuts
1095 if (fUseMCInfo) {
1096
1097 // find associated MC particle for Lc -> p+K+pi
1098 Int_t mcLabel = d->MatchToMC(4122, mcArray, 3, pdgDgLctopKpi);
1099 //Int_t mcLabelTemp = d->MatchToMC(4122, mcArray);
1100 //Printf("mcLabel = %d, mcLabelTemp = %d", mcLabel, mcLabelTemp);
1101 if (mcLabel >= 0) {
1102
1103 AliAODMCParticle *partLcpKpi = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcLabel));
1104 if(partLcpKpi){
1105 Int_t pdgCode = partLcpKpi->GetPdgCode();
1106 AliDebug(2,Form(" ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ MClabel=%d ~~~~~~~~~~ pdgCode=%d", mcLabel, pdgCode));
1107 //AliInfo(Form(" ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ MClabel=%d ~~~~~~~~~~ pdgCode=%d", mcLabel, pdgCode));
1108 fHistoLcpKpiBeforeCuts->Fill(1);
1109
1110 }
1111 }
1112 else {
1113 //AliInfo(Form(" ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ MClabel=%d ~~~~~~~~~~", mcLabel));
1114 fHistoLcpKpiBeforeCuts->Fill(0);
1115 }
1116 }
1117 }
1118
1119 // loop over cascades to search for candidates Lc->p+K0S
1120
1121 for (Int_t iLctopK0s = 0; iLctopK0s < nCascades; iLctopK0s++) {
1122
1123 // Lc candidates and K0s from Lc
1124 AliAODRecoCascadeHF* lcK0spr = dynamic_cast<AliAODRecoCascadeHF*>(arrayLctopKos->At(iLctopK0s));
1125 if (!lcK0spr) {
1126 AliDebug(2,Form("Cascade %d doens't exist, skipping",iLctopK0s));
1127 continue;
1128 }
1129
1130 //Filling a control histogram with no cuts
1131 if (fUseMCInfo) {
1132
1133 Int_t pdgCode=-2;
1134
1135 // find associated MC particle for Lc -> p+K0 and K0S->pi+pi
1136 fmcLabelLc = lcK0spr->MatchToMC(pdgCand, pdgDgLctoV0bachelor[1], pdgDgLctoV0bachelor, pdgDgV0toDaughters, mcArray, kTRUE);
1137 if (fmcLabelLc>=0) {
4539d15b 1138 AliDebug(2, Form("----> cascade number %d (total cascade number = %d) is a Lc!", iLctopK0s, nCascades));
0abe5247 1139
1140 AliAODMCParticle *partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(fmcLabelLc));
1141 if(partLc){
1142 pdgCode = partLc->GetPdgCode();
1143 if (pdgCode<0) AliDebug(2,Form(" ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ MClabel=%d ~~~~~~~~~~ pdgCode=%d", fmcLabelLc, pdgCode));
1144 pdgCode = TMath::Abs(pdgCode);
1145 fHistoLcBeforeCuts->Fill(1);
1146
1147 }
1148 }
1149 else {
1150 fHistoLcBeforeCuts->Fill(0);
1151 pdgCode=-1;
1152 }
1153 }
1154
1155 //if (!lcK0spr->GetSecondaryVtx()) {
1156 // AliInfo("No secondary vertex");
1157 //continue;
1158 //}
1159
1160 if (lcK0spr->GetNDaughters()!=2) {
1161 AliDebug(2,Form("Cascade %d has not 2 daughters (nDaughters=%d)",iLctopK0s,lcK0spr->GetNDaughters()));
1162 continue;
1163 }
1164
1165 AliAODv0 * v0part = dynamic_cast<AliAODv0*>(lcK0spr->Getv0());
1166 AliAODTrack * bachPart = dynamic_cast<AliAODTrack*>(lcK0spr->GetBachelor());
1167 if (!v0part || !bachPart) {
1168 AliDebug(2,Form("Cascade %d has no V0 or no bachelor object",iLctopK0s));
1169 continue;
1170 }
1171
1172
1173 if (!v0part->GetSecondaryVtx()) {
1174 AliDebug(2,Form("No secondary vertex for V0 by cascade %d",iLctopK0s));
1175 continue;
1176 }
1177
1178 if (v0part->GetNDaughters()!=2) {
1179 AliDebug(2,Form("current V0 has not 2 daughters (onTheFly=%d, nDaughters=%d)",v0part->GetOnFlyStatus(),v0part->GetNDaughters()));
1180 continue;
1181 }
1182
1183 AliAODTrack * v0Pos = dynamic_cast<AliAODTrack*>(lcK0spr->Getv0PositiveTrack());
1184 AliAODTrack * v0Neg = dynamic_cast<AliAODTrack*>(lcK0spr->Getv0NegativeTrack());
ff22c824 1185 if (!v0Neg || !v0Pos) {
0abe5247 1186 AliDebug(2,Form("V0 by cascade %d has no V0positive of V0negative object",iLctopK0s));
1187 continue;
1188 }
1189
1190
1191 if (v0Pos->Charge() == v0Neg->Charge()) {
1192 AliDebug(2,Form("V0 by cascade %d has daughters with the same sign: IMPOSSIBLE!",iLctopK0s));
1193 continue;
1194 }
1195
1196 Int_t isLc = 0;
1197
1198 if (fUseMCInfo) {
1199
1200 Int_t pdgCode=-2;
1201
1202 // find associated MC particle for Lc -> p+K0 and K0S->pi+pi
1203 Int_t mcLabel = lcK0spr->MatchToMC(pdgCand, pdgDgLctoV0bachelor[1], pdgDgLctoV0bachelor, pdgDgV0toDaughters, mcArray, kTRUE);
1204 if (mcLabel>=0) {
1205 AliDebug(2,Form(" ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~cascade number %d (total cascade number = %d)", iLctopK0s, nCascades));
1206
1207 AliAODMCParticle *partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcLabel));
1208 if(partLc){
1209 pdgCode = partLc->GetPdgCode();
1210 if (pdgCode<0) AliDebug(2,Form(" ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ MClabel=%d ~~~~~~~~~~ pdgCode=%d", mcLabel, pdgCode));
1211 pdgCode = TMath::Abs(pdgCode);
1212 isLc = 1;
1213 fHistoLc->Fill(1);
1214
1215 }
1216 }
1217 else {
1218 fHistoLc->Fill(0);
1219 pdgCode=-1;
1220 }
1221 }
1222 AliDebug(2, Form("\n\n\n Analysing candidate %d\n", iLctopK0s));
1223 AliDebug(2, Form(">>>>>>>>>> Candidate is background, fFillOnlySgn = %d --> SKIPPING", fFillOnlySgn));
1224 if (!isLc) {
1225 if (fFillOnlySgn) { // if it is background, and we want only signal, we do not fill the tree
1226 continue;
1227 }
1228 else { // checking if we want to fill the background
1229 if (fKeepingOnlyHIJINGBkg){
1230 // we have decided to fill the background only when the candidate has the daugthers that all come from HIJING underlying event!
122f753d 1231 Bool_t isCandidateInjected = fUtils->HasCascadeCandidateAnyDaughInjected(lcK0spr, aodheader, mcArray);
0abe5247 1232 if (!isCandidateInjected){
1233 AliDebug(2, "The candidate is from HIJING (i.e. not injected), keeping it to fill background");
1234 fHistoBackground->Fill(1);
1235 }
1236 else {
1237 AliDebug(2, "The candidate is NOT from HIJING, we skip it when filling background");
1238 fHistoBackground->Fill(0);
1239 continue;
1240 }
1241 }
ae48be26 1242 else if (fKeepingOnlyPYTHIABkg){
1243 // we have decided to fill the background only when the candidate has the daugthers that all come from HIJING underlying event!
1244 AliAODTrack *bachelor = (AliAODTrack*)lcK0spr->GetBachelor();
1245 AliAODTrack *v0pos = (AliAODTrack*)lcK0spr->Getv0PositiveTrack();
1246 AliAODTrack *v0neg = (AliAODTrack*)lcK0spr->Getv0NegativeTrack();
1247 if (!bachelor || !v0pos || !v0neg) {
1248 AliDebug(2, "Cannot retrieve one of the tracks while checking origin, continuing");
1249 continue;
1250 }
1251 else {
1252 Int_t labelbachelor = TMath::Abs(bachelor->GetLabel());
1253 Int_t labelv0pos = TMath::Abs(v0pos->GetLabel());
1254 Int_t labelv0neg = TMath::Abs(v0neg->GetLabel());
1255 AliAODMCParticle* MCbachelor = (AliAODMCParticle*)mcArray->At(labelbachelor);
1256 AliAODMCParticle* MCv0pos = (AliAODMCParticle*)mcArray->At(labelv0pos);
1257 AliAODMCParticle* MCv0neg = (AliAODMCParticle*)mcArray->At(labelv0neg);
1258 if (!MCbachelor || !MCv0pos || !MCv0neg) {
1259 AliDebug(2, "Cannot retrieve MC particle for one of the tracks while checking origin, continuing");
1260 continue;
1261 }
1262 else {
1263 Int_t isBachelorFromPythia = fUtils->CheckOrigin(mcArray, MCbachelor, kTRUE);
1264 Int_t isv0posFromPythia = fUtils->CheckOrigin(mcArray, MCv0pos, kTRUE);
1265 Int_t isv0negFromPythia = fUtils->CheckOrigin(mcArray, MCv0neg, kTRUE);
1266 if (isBachelorFromPythia != 0 && isv0posFromPythia != 0 && isv0negFromPythia != 0){
1267 AliDebug(2, "The candidate is from PYTHIA (i.e. all daughters originate from a quark), keeping it to fill background");
1268 fHistoBackground->Fill(2);
1269 }
1270 else {
1271 AliDebug(2, "The candidate is NOT from PYTHIA, we skip it when filling background");
1272 fHistoBackground->Fill(3);
1273 continue;
1274 }
1275 }
1276 }
1277 }
0abe5247 1278 }
1279 }
1280
44c36327 1281 FillLc2pK0Sspectrum(lcK0spr, isLc, nSelectedAnal, cutsAnal, mcArray, iLctopK0s);
0abe5247 1282
1283 }
1284
1285 return;
1286
1287}
1288//________________________________________________________________________
1289void AliAnalysisTaskSELc2V0bachelorTMVA::FillLc2pK0Sspectrum(AliAODRecoCascadeHF *part,
1290 Int_t isLc,
1291 Int_t &nSelectedAnal,
1292 AliRDHFCutsLctoV0 *cutsAnal,
44c36327 1293 TClonesArray *mcArray, Int_t iLctopK0s){
0abe5247 1294 //
1295 // Fill histos for Lc -> K0S+proton
1296 //
1297
1298 /*
1299 if (!part->GetOwnPrimaryVtx()) {
1300 //Printf("No primary vertex for Lc found!!");
1301 part->SetOwnPrimaryVtx(fVtx1);
1302 }
1303 else {
1304 //Printf("Yu-huuuu!!! primary vertex for Lc found!!");
1305 }
1306 */
1307 Double_t invmassLc = part->InvMassLctoK0sP();
1308
1309 AliAODv0 * v0part = part->Getv0();
1310 Bool_t onFlyV0 = v0part->GetOnFlyStatus(); // on-the-flight V0s
1311 if (onFlyV0){ // on-the-fly V0
1312 if (isLc) { // Lc
1313 fHistoLcOnTheFly->Fill(2.);
1314 }
1315 else { // not Lc
1316 fHistoLcOnTheFly->Fill(0.);
1317 }
1318 }
1319 else { // offline V0
1320 if (isLc) { // Lc
1321 fHistoLcOnTheFly->Fill(3.);
1322 }
1323 else { // not Lc
1324 fHistoLcOnTheFly->Fill(1.);
1325 }
1326 }
1327
1328 Double_t dcaV0 = v0part->GetDCA();
1329 Double_t invmassK0s = v0part->MassK0Short();
1330
1331 if ( (cutsAnal->IsInFiducialAcceptance(part->Pt(),part->Y(4122))) ) {
1332 if (isLc) {
1333 fHistoFiducialAcceptance->Fill(3.);
1334 }
1335 else {
1336 fHistoFiducialAcceptance->Fill(1.);
1337 }
1338 }
1339 else {
1340 if (isLc) {
1341 fHistoFiducialAcceptance->Fill(2.);
1342 }
1343 else {
1344 fHistoFiducialAcceptance->Fill(0.);
1345 }
1346 }
1347
1348 Int_t isInV0window = (((cutsAnal->IsSelectedSingleCut(part, AliRDHFCuts::kCandidate, 2)) & (AliRDHFCutsLctoV0::kLcToK0Spr)) == (AliRDHFCutsLctoV0::kLcToK0Spr)); // cut on V0 invMass
1349
1350 if (isInV0window == 0) {
44c36327 1351 AliDebug(2, "No: The candidate has NOT passed the V0 window cuts!");
1352 if (isLc) Printf("SIGNAL candidate rejected: V0 window cuts");
0abe5247 1353 return;
1354 }
1355 else AliDebug(2, "Yes: The candidate has passed the mass cuts!");
1356
0abe5247 1357 Bool_t isInCascadeWindow = (((cutsAnal->IsSelectedSingleCut(part, AliRDHFCuts::kCandidate, 0)) & (AliRDHFCutsLctoV0::kLcToK0Spr)) == (AliRDHFCutsLctoV0::kLcToK0Spr)); // cut on Lc->p+K0S invMass
1358
1359 if (!isInCascadeWindow) {
1360 AliDebug(2, "No: The candidate has NOT passed the cascade window cuts!");
44c36327 1361 if (isLc) Printf("SIGNAL candidate rejected: cascade window cuts");
0abe5247 1362 return;
1363 }
1364 else AliDebug(2, "Yes: The candidate has passed the cascade window cuts!");
1365
1366 Bool_t isCandidateSelectedCuts = (((cutsAnal->IsSelected(part, AliRDHFCuts::kCandidate)) & (AliRDHFCutsLctoV0::kLcToK0Spr)) == (AliRDHFCutsLctoV0::kLcToK0Spr)); // kinematic/topological cuts
44c36327 1367 AliDebug(2, Form("recoAnalysisCuts = %d", cutsAnal->IsSelected(part, AliRDHFCuts::kCandidate) & (AliRDHFCutsLctoV0::kLcToK0Spr)));
0abe5247 1368 if (!isCandidateSelectedCuts){
1369 AliDebug(2, "No: Analysis cuts kCandidate level NOT passed");
1370 if (isLc) Printf("SIGNAL candidate rejected");
1371 return;
1372 }
1373 else {
1374 AliDebug(2, "Yes: Analysis cuts kCandidate level passed");
1375 }
1376
1377 AliAODTrack *bachelor = (AliAODTrack*)part->GetBachelor();
1378 if (!bachelor) {
1379 AliDebug(2, Form("Very weird, the bachelor is not there... returning for this candidate"));
1380 return;
1381 }
1382
1383 //Bool_t isBachelorID = (((cutsAnal->IsSelected(part,AliRDHFCuts::kPID))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)); // ID x bachelor
ae48be26 1384 Double_t probTPCTOF[AliPID::kSPECIES]={-1.};
0abe5247 1385
1386 UInt_t detUsed = fPIDCombined->ComputeProbabilities(bachelor, fPIDResponse, probTPCTOF);
2b0b36b0 1387 AliDebug(2, Form("detUsed (TPCTOF case) = %d", detUsed));
ae48be26 1388 Double_t probProton = -1.;
6be55058 1389 // Double_t probPion = -1.;
1390 // Double_t probKaon = -1.;
0abe5247 1391 if (detUsed == (UInt_t)fPIDCombined->GetDetectorMask() ) {
2b0b36b0 1392 AliDebug(2, Form("We have found the detector mask for TOF + TPC: probProton will be set to %f", probTPCTOF[AliPID::kProton]));
0abe5247 1393 probProton = probTPCTOF[AliPID::kProton];
6be55058 1394 // probPion = probTPCTOF[AliPID::kPion];
1395 // probKaon = probTPCTOF[AliPID::kKaon];
0abe5247 1396 }
ae48be26 1397 else { // if you don't have both TOF and TPC, try only TPC
1398 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC);
2b0b36b0 1399 AliDebug(2, "We did not find the detector mask for TOF + TPC, let's see only TPC");
ae48be26 1400 detUsed = fPIDCombined->ComputeProbabilities(bachelor, fPIDResponse, probTPCTOF);
2b0b36b0 1401 AliDebug(2,Form(" detUsed (TPC case) = %d", detUsed));
ae48be26 1402 if (detUsed == (UInt_t)fPIDCombined->GetDetectorMask()) {
1403 probProton = probTPCTOF[AliPID::kProton];
6be55058 1404 // probPion = probTPCTOF[AliPID::kPion];
1405 // probKaon = probTPCTOF[AliPID::kKaon];
2b0b36b0 1406 AliDebug(2, Form("TPC only worked: probProton will be set to %f", probTPCTOF[AliPID::kProton]));
ae48be26 1407 }
1408 else {
2b0b36b0 1409 AliDebug(2, "Only TPC did not work...");
ae48be26 1410 }
1411 // resetting mask to ask for both TPC+TOF
1412 fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC+AliPIDResponse::kDetTOF);
1413 }
2b0b36b0 1414 AliDebug(2, Form("probProton = %f", probProton));
ae48be26 1415
1416 // now we get the TPC and TOF single PID probabilities (only for Proton, or the tree will explode :) )
1417 Double_t probProtonTPC = -1.;
1418 Double_t probProtonTOF = -1.;
1419 Double_t pidTPC[AliPID::kSPECIES]={-1.};
1420 Double_t pidTOF[AliPID::kSPECIES]={-1.};
1421 Int_t respTPC = fPIDResponse->ComputePIDProbability(AliPIDResponse::kDetTPC, bachelor, AliPID::kSPECIES, pidTPC);
1422 Int_t respTOF = fPIDResponse->ComputePIDProbability(AliPIDResponse::kDetTOF, bachelor, AliPID::kSPECIES, pidTOF);
1423 if (respTPC == AliPIDResponse::kDetPidOk) probProtonTPC = pidTPC[AliPID::kProton];
1424 if (respTOF == AliPIDResponse::kDetPidOk) probProtonTOF = pidTOF[AliPID::kProton];
1425
1426 // checking V0 status (on-the-fly vs offline)
0abe5247 1427 if ( !( !onFlyV0 || (onFlyV0 && fUseOnTheFlyV0) ) ) {
1428 AliDebug(2, "On-the-fly discarded");
1429 return;
1430 }
1431
1432 if ( (((cutsAnal->IsSelected(part,AliRDHFCuts::kAll))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)) ) {
1433 nSelectedAnal++;
1434 }
1435
1436 if ( !(cutsAnal->IsInFiducialAcceptance(part->Pt(),part->Y(4122))) ) return;
1437
1438 if ( !( ( (cutsAnal->IsSelected(part, AliRDHFCuts::kTracks)) & (AliRDHFCutsLctoV0::kLcToK0Spr)) == (AliRDHFCutsLctoV0::kLcToK0Spr) ) ) { // esd track cuts
1439 if (isLc) Printf("SIGNAL candidate rejected");
1440 AliDebug(2, "No: Analysis cuts kTracks level NOT passed");
1441 return;
1442 }
1443 else {
1444 AliDebug(2, "Yes: Analysis cuts kTracks level passed");
1445 }
1446
1447 Int_t pdgCand = 4122;
1448 Int_t pdgDgLctoV0bachelor[2]={211, 3122}; // case of wrong decay! Lc --> L + pi
1449 Int_t pdgDgV0toDaughters[2]={2212, 211}; // case of wrong decay! Lc --> L + pi
1450 Int_t isLc2LBarpi=0, isLc2Lpi=0;
1451 Int_t currentLabel = part->GetLabel();
1452 Int_t mcLabel = 0;
1453 if (fUseMCInfo) {
1454 mcLabel = part->MatchToMC(pdgCand, pdgDgLctoV0bachelor[1], pdgDgLctoV0bachelor, pdgDgV0toDaughters, mcArray, kTRUE);
1455 if (mcLabel>=0) {
1456 if (bachelor->Charge()==-1) isLc2LBarpi=1;
1457 if (bachelor->Charge()==+1) isLc2Lpi=1;
1458 }
1459 }
1460
1461 Int_t pdgDg2prong[2] = {211, 211};
1462 Int_t labelK0S = 0;
1463 Int_t isK0S = 0;
1464 if (fUseMCInfo) {
1465 labelK0S = v0part->MatchToMC(310, mcArray, 2, pdgDg2prong);
1466 if (labelK0S>=0) isK0S = 1;
1467 }
1468
1469 pdgDg2prong[0] = 211;
1470 pdgDg2prong[1] = 2212;
1471 Int_t isLambda = 0;
1472 Int_t isLambdaBar = 0;
1473 Int_t lambdaLabel = 0;
1474 if (fUseMCInfo) {
1475 lambdaLabel = v0part->MatchToMC(3122, mcArray, 2, pdgDg2prong);
1476 if (lambdaLabel>=0) {
1477 AliAODMCParticle *lambdaTrack = (AliAODMCParticle*)mcArray->At(lambdaLabel);
1478 if (lambdaTrack->GetPdgCode()==3122) isLambda = 1;
1479 else if (lambdaTrack->GetPdgCode()==-3122) isLambdaBar = 1;
1480 }
1481 }
1482
1483 pdgDg2prong[0] = 11;
1484 pdgDg2prong[1] = 11;
1485 Int_t isGamma = 0;
1486 Int_t gammaLabel = 0;
1487 if (fUseMCInfo) {
1488 gammaLabel = v0part->MatchToMC(22, mcArray, 2, pdgDg2prong);
1489 if (gammaLabel>=0) {
1490 AliAODMCParticle *gammaTrack = (AliAODMCParticle*)mcArray->At(gammaLabel);
1491 if (gammaTrack->GetPdgCode()==22) isGamma = 1;
1492 }
1493 }
1494
1495 Int_t pdgTemp = -1;
1496 if (currentLabel != -1){
1497 AliAODMCParticle *tempPart = (AliAODMCParticle*)mcArray->At(currentLabel);
1498 pdgTemp = tempPart->GetPdgCode();
1499 }
1500 if (isLc) AliDebug(2, Form("Signal: Candidate is a Lc in K0s+p"));
1501 else if (isLc2LBarpi) AliDebug(2, Form("Background: Candidate is a Lc in Lbar + pi"));
1502 else if (isLc2Lpi) AliDebug(2, Form("Background: Candidate is a Lc in L + pi"));
1503 else AliDebug(2, Form("Pure bkg: Candidate has pdg = %d", pdgTemp));
1504 if (isK0S) AliDebug(2, Form("V0 is a K0S"));
1505 else if (isLambda) AliDebug(2, Form("V0 is a Lambda"));
1506 else if (isLambdaBar) AliDebug(2, Form("V0 is a LambdaBar"));
1507 else if (isGamma) AliDebug(2, Form("V0 is a Gamma"));
1508 //else AliDebug(2, Form("V0 is something else!!"));
1509
1510 Double_t invmassLc2Lpi = part->InvMassLctoLambdaPi();
1511 Double_t invmassLambda = v0part->MassLambda();
1512 Double_t invmassLambdaBar = v0part->MassAntiLambda();
1513
1514 //Double_t nSigmaITSpr=-999.;
1515 Double_t nSigmaTPCpr=-999.;
1516 Double_t nSigmaTOFpr=-999.;
1517
1518 //Double_t nSigmaITSpi=-999.;
1519 Double_t nSigmaTPCpi=-999.;
1520 Double_t nSigmaTOFpi=-999.;
1521
1522 //Double_t nSigmaITSka=-999.;
1523 Double_t nSigmaTPCka=-999.;
1524 Double_t nSigmaTOFka=-999.;
1525
1526 /*
1527 cutsAnal->GetPidHF()->GetnSigmaITS(bachelor,4,nSigmaITSpr);
1528 cutsAnal->GetPidHF()->GetnSigmaTPC(bachelor,4,nSigmaTPCpr);
1529 cutsAnal->GetPidHF()->GetnSigmaTOF(bachelor,4,nSigmaTOFpr);
1530 cutsAnal->GetPidHF()->GetnSigmaITS(bachelor,2,nSigmaITSpi);
1531 cutsAnal->GetPidHF()->GetnSigmaTPC(bachelor,2,nSigmaTPCpi);
1532 cutsAnal->GetPidHF()->GetnSigmaTOF(bachelor,2,nSigmaTOFpi);
1533 cutsAnal->GetPidHF()->GetnSigmaITS(bachelor,3,nSigmaITSka);
1534 cutsAnal->GetPidHF()->GetnSigmaTPC(bachelor,3,nSigmaTPCka);
1535 cutsAnal->GetPidHF()->GetnSigmaTOF(bachelor,3,nSigmaTOFka);
1536 */
1537
1538 nSigmaTPCpi = fPIDResponse->NumberOfSigmasTPC(bachelor,(AliPID::kPion));
1539 nSigmaTPCka = fPIDResponse->NumberOfSigmasTPC(bachelor,(AliPID::kKaon));
1540 nSigmaTPCpr = fPIDResponse->NumberOfSigmasTPC(bachelor,(AliPID::kProton));
1541 nSigmaTOFpi = fPIDResponse->NumberOfSigmasTOF(bachelor,(AliPID::kPion));
1542 nSigmaTOFka = fPIDResponse->NumberOfSigmasTOF(bachelor,(AliPID::kKaon));
1543 nSigmaTOFpr = fPIDResponse->NumberOfSigmasTOF(bachelor,(AliPID::kProton));
1544
1545
1546 // Fill candidate variable Tree (track selection, V0 invMass selection)
1547 if (!onFlyV0 && isInV0window && isInCascadeWindow && part->CosV0PointingAngle()>0.99 && TMath::Abs(nSigmaTPCpr) <= 3 && v0part->Getd0Prong(0) < 20 && v0part->Getd0Prong(1) < 20) {
1548
1549 fCandidateVariables[0] = invmassLc;
1550 fCandidateVariables[1] = invmassLc2Lpi;
1551 fCandidateVariables[2] = invmassK0s;
1552 fCandidateVariables[3] = invmassLambda;
1553 fCandidateVariables[4] = invmassLambdaBar;
1554 fCandidateVariables[5] = part->CosV0PointingAngle();
1555 fCandidateVariables[6] = dcaV0;
1556 fCandidateVariables[7] = part->Getd0Prong(0);
1557 fCandidateVariables[8] = part->Getd0Prong(1);
1558 fCandidateVariables[9] = nSigmaTPCpr;
1559 fCandidateVariables[10] = nSigmaTPCpi;
1560 fCandidateVariables[11] = nSigmaTPCka;
1561 fCandidateVariables[12] = nSigmaTOFpr;
1562 fCandidateVariables[13] = nSigmaTOFpi;
1563 fCandidateVariables[14] = nSigmaTOFka;
1564 fCandidateVariables[15] = bachelor->Pt();
1565 AliAODTrack *v0pos = (AliAODTrack*)part->Getv0PositiveTrack();
1566 fCandidateVariables[16] = v0pos->Pt();
1567 AliAODTrack *v0neg = (AliAODTrack*)part->Getv0NegativeTrack();
1568 fCandidateVariables[17] = v0neg->Pt();
1569 fCandidateVariables[18] = v0part->Getd0Prong(0);
1570 fCandidateVariables[19] = v0part->Getd0Prong(1);
1571 fCandidateVariables[20] = v0part->Pt();
1572 fCandidateVariables[21] = v0part->InvMass2Prongs(0,1,11,11);
1573 fCandidateVariables[22] = part->Pt();
1574 fCandidateVariables[23] = probProton;
1575 fCandidateVariables[24] = part->Eta();
1576 fCandidateVariables[25] = v0pos->Eta();
1577 fCandidateVariables[26] = v0neg->Eta();
ae48be26 1578 fCandidateVariables[27] = probProtonTPC;
1579 fCandidateVariables[28] = probProtonTOF;
0abe5247 1580 fCandidateVariables[29] = bachelor->Eta();
1581
1582 fCandidateVariables[30] = part->P();
1583 fCandidateVariables[31] = bachelor->P();
1584 fCandidateVariables[32] = v0part->P();
1585 fCandidateVariables[33] = v0pos->P();
1586 fCandidateVariables[34] = v0neg->P();
1587
1588 fCandidateVariables[35] = part->Y(4122);
1589 fCandidateVariables[36] = bachelor->Y(2212);
1590 fCandidateVariables[37] = v0part->Y(310);
1591 fCandidateVariables[38] = v0pos->Y(211);
1592 fCandidateVariables[39] = v0neg->Y(211);
1593
1594 fCandidateVariables[40] = v0part->Eta();
1595
1596 fCandidateVariables[41] = part->DecayLength();
1597 fCandidateVariables[42] = part->DecayLengthV0();
1598 fCandidateVariables[43] = part->Ct(4122);
1599 fCandidateVariables[44] = v0part->Ct(310, v0part->GetSecondaryVtx());
1600
ef87f0d6 1601 EBachelor bachCode = kBachInvalid;
1602 EK0S k0SCode = kK0SInvalid;
1603 if (fUseMCInfo) {
1604 bachCode = CheckBachelor(part, bachelor, mcArray);
1605 k0SCode = CheckK0S(part, v0part, mcArray);
1606 }
0abe5247 1607
1608 fCandidateVariables[45] = bachCode;
1609 fCandidateVariables[46] = k0SCode;
1610
1611 Double_t V0KF[3] = {-999999, -999999, -999999};
1612 Double_t errV0KF[3] = {-999999, -999999, -999999};
1613 Double_t LcKF[3] = {-999999, -999999, -999999};
1614 Double_t errLcKF[3] = {-999999, -999999, -999999};
1615 Double_t distances[3] = {-999999, -999999, -999999};
1616 Double_t armPolKF[2] = {-999999, -999999};
1617
1618 if (fCallKFVertexing){
1619 Int_t kfResult = CallKFVertexing(part, v0part, bachelor, mcArray, &V0KF[0], &errV0KF[0], &LcKF[0], &errLcKF[0], &distances[0], &armPolKF[0]);
1620 AliDebug(2, Form("Result from KF = %d", kfResult));
1621 }
1622
1623 /*
1624 for (Int_t i = 0; i< 3; i++){
1625 Printf("i = %d, V0KF = %f, errV0KF = %f, LcKF = %f, errLcKF = %f", V0KF[i], errV0KF[i], LcKF[i], errLcKF[i]);
1626 }
1627 */
1628
1629 fCandidateVariables[47] = V0KF[0];
1630 fCandidateVariables[48] = V0KF[1];
1631 fCandidateVariables[49] = V0KF[2];
1632
1633 fCandidateVariables[50] = errV0KF[0];
1634 fCandidateVariables[51] = errV0KF[1];
1635 fCandidateVariables[52] = errV0KF[2];
1636
1637 fCandidateVariables[53] = LcKF[0];
1638 fCandidateVariables[54] = LcKF[1];
1639 fCandidateVariables[55] = LcKF[2];
1640
1641 fCandidateVariables[56] = errLcKF[0];
1642 fCandidateVariables[57] = errLcKF[1];
1643 fCandidateVariables[58] = errLcKF[2];
1644
1645 fCandidateVariables[59] = distances[0];
1646 fCandidateVariables[60] = distances[1];
1647 fCandidateVariables[61] = distances[2];
1648 fCandidateVariables[62] = armPolKF[0];
1649 fCandidateVariables[63] = armPolKF[1];
1650 fCandidateVariables[64] = v0part->AlphaV0();
1651 fCandidateVariables[65] = v0part->PtArmV0();
1652
ef87f0d6 1653 AliDebug(2, Form("v0pos->GetStatus() & AliESDtrack::kITSrefit= %d, v0neg->GetStatus() & AliESDtrack::kITSrefit = %d, v0pos->GetTPCClusterInfo(2, 1)= %f, v0neg->GetTPCClusterInfo(2, 1) = %f", (Int_t)(v0pos->GetStatus() & AliESDtrack::kITSrefit), (Int_t)(v0pos->GetStatus() & AliESDtrack::kITSrefit), v0pos->GetTPCClusterInfo(2, 1), v0neg->GetTPCClusterInfo(2, 1)));
1654 fCandidateVariables[66] = v0pos->GetStatus() & AliESDtrack::kITSrefit;
1655 fCandidateVariables[67] = v0neg->GetStatus() & AliESDtrack::kITSrefit;
1656 fCandidateVariables[68] = v0pos->GetTPCClusterInfo(2, 1);
1657 fCandidateVariables[69] = v0neg->GetTPCClusterInfo(2, 1);
1658
1659 fCandidateVariables[70] = v0part->Xv();
1660 fCandidateVariables[71] = v0part->Yv();
1661 fCandidateVariables[72] = v0part->Zv();
1662
1663 fCandidateVariables[73] = fVtx1->GetX();
1664 fCandidateVariables[74] = fVtx1->GetY();
1665 fCandidateVariables[75] = fVtx1->GetZ();
1666
1667 fCandidateVariables[76] = bachelor->GetITSNcls();
1668 fCandidateVariables[77] = bachelor->HasPointOnITSLayer(0) + bachelor->HasPointOnITSLayer(1);
1669
1670 fCandidateVariables[78] = v0pos->GetITSNcls();
1671 fCandidateVariables[79] = v0pos->HasPointOnITSLayer(0) + v0pos->HasPointOnITSLayer(1);
1672
1673 fCandidateVariables[80] = v0neg->GetITSNcls();
1674 fCandidateVariables[81] = v0neg->HasPointOnITSLayer(0) + v0neg->HasPointOnITSLayer(1);
1675
1676 TVector3 mom1(bachelor->Px(), bachelor->Py(), bachelor->Pz());
1677 TVector3 mom2(v0part->Px(), v0part->Py(), v0part->Pz());
1678 TVector3 momTot(part->Px(), part->Py(), part->Pz());
1679
1680 Double_t Ql1 = mom1.Dot(momTot)/momTot.Mag();
1681 Double_t Ql2 = mom2.Dot(momTot)/momTot.Mag();
1682
1683 Double_t alphaArmLc = (Ql1 - Ql2)/(Ql1 + Ql2);
1684 Double_t alphaArmLcCharge = ( bachelor->Charge() > 0 ? (Ql1 - Ql2)/(Ql1 + Ql2) : (Ql2 - Ql1)/(Ql1 + Ql2) );
1685 Double_t ptArmLc = mom1.Perp(momTot);
1686
1687 fCandidateVariables[82] = alphaArmLc;
1688 fCandidateVariables[83] = alphaArmLcCharge;
1689 fCandidateVariables[84] = ptArmLc;
1690
1691 Double_t massK0SPDG = TDatabasePDG::Instance()->GetParticle(310)->Mass(); // mass K0S PDG
1692 Double_t massPrPDG = TDatabasePDG::Instance()->GetParticle(2212)->Mass(); // mass Proton PDG
1693 Double_t massLcPDG = TDatabasePDG::Instance()->GetParticle(4122)->Mass(); // mass Lc PDG
1694
1695 Double_t pStar = TMath::Sqrt((massLcPDG*massLcPDG-massPrPDG*massPrPDG-massK0SPDG*massK0SPDG)*(massLcPDG*massLcPDG-massPrPDG*massPrPDG-massK0SPDG*massK0SPDG)-4.*massPrPDG*massPrPDG*massK0SPDG*massK0SPDG)/(2.*massLcPDG);
1696 Double_t e = part->E(4122);
1697 Double_t beta = part->P()/e;
1698 Double_t gamma = e/massLcPDG;
1699
1700 Double_t cts = (Ql1/gamma-beta*TMath::Sqrt(pStar*pStar+massPrPDG*massPrPDG))/pStar;
1701
1702 fCandidateVariables[85] = cts;
1703
0abe5247 1704 if (fUseMCInfo) {
1705 if (isLc){
44c36327 1706 AliDebug(2, Form("Reco particle %d --> Filling Sgn", iLctopK0s));
0abe5247 1707 fVariablesTreeSgn->Fill();
1708 fHistoCodesSgn->Fill(bachCode, k0SCode);
1709 }
1710 else {
1711 if (fFillOnlySgn == kFALSE){
1712 AliDebug(2, "Filling Bkg");
1713 fVariablesTreeBkg->Fill();
1714 fHistoCodesBkg->Fill(bachCode, k0SCode);
1715 }
1716 }
1717 }
1718 else {
1719 fVariablesTreeSgn->Fill();
1720 }
1721 }
1722
1723 return;
1724
1725}
1726
1727//________________________________________________________________________
1728Int_t AliAnalysisTaskSELc2V0bachelorTMVA::CallKFVertexing(AliAODRecoCascadeHF *cascade, AliAODv0* v0part, AliAODTrack* bach, TClonesArray *mcArray,
1729 Double_t* V0KF, Double_t* errV0KF, Double_t* LcKF, Double_t* errLcKF,
1730 Double_t* distances, Double_t* armPolKF) {
1731
1732 //
1733 // method to perform KF vertexing
1734 // elements: [0] = mass, [1] = DecayLength, [2] = lifeTime
1735 //
1736
1737 Int_t codeKFV0 = -1, codeKFLc = -1;
1738
1739 AliKFVertex primVtxCopy;
1740 Int_t nt = 0, ntcheck = 0;
1741 Double_t pos[3] = {0., 0., 0.};
1742
1743 fVtx1->GetXYZ(pos);
1744 Int_t contr = fVtx1->GetNContributors();
1745 Double_t covmatrix[6] = {0.};
1746 fVtx1->GetCovarianceMatrix(covmatrix);
1747 Double_t chi2 = fVtx1->GetChi2();
1748 AliESDVertex primaryESDVtxCopy(pos, covmatrix, chi2, contr, "Vertex");
1749
1750 // topological constraint
1751 primVtxCopy = AliKFVertex(primaryESDVtxCopy);
1752 nt = primaryESDVtxCopy.GetNContributors();
1753 ntcheck = nt;
1754
1755 Int_t pdg[2] = {211, -211};
1756 Int_t pdgLc[2] = {2212, 310};
1757
1758 Int_t pdgDgV0toDaughters[2] = {211, 211};
1759
1760 Int_t mcLabelV0 = v0part->MatchToMC(310, mcArray, 2, pdgDgV0toDaughters);
1761
1762 // the KF vertex for the V0 has to be built with the prongs of the V0!
1763 Bool_t isMCokV0 = kTRUE, isBkgV0 = kFALSE;
ef87f0d6 1764 AliKFParticle V0, positiveV0KF, negativeV0KF;
0abe5247 1765 Int_t labelsv0daugh[2] = {-1, -1};
ef87f0d6 1766 Int_t idv0daugh[2] = {-1, -1};
1767 AliExternalTrackParam* esdv0Daugh1 = 0x0;
1768 AliExternalTrackParam* esdv0Daugh2 = 0x0;
0abe5247 1769 for(Int_t ipr= 0; ipr < 2; ipr++){ // 0 is positive, 1 is negative
1770 AliAODTrack *aodTrack = (AliAODTrack*)v0part->GetDaughter(ipr);
1771 if(!aodTrack) {
1772 AliDebug(2, "No V0 daughters available");
1773 return -1;
1774 }
ef87f0d6 1775 Double_t xyz[3], pxpypz[3], cv[21];
1776 Short_t sign;
1777 aodTrack->GetXYZ(xyz);
1778 aodTrack->PxPyPz(pxpypz);
1779 aodTrack->GetCovarianceXYZPxPyPz(cv);
1780 sign = aodTrack->Charge();
1781 AliExternalTrackParam tmp1( xyz, pxpypz, cv, sign);
1782
1783 if (ipr == 0) esdv0Daugh1 = new AliExternalTrackParam( xyz, pxpypz, cv, sign);
1784 else esdv0Daugh2 = new AliExternalTrackParam( xyz, pxpypz, cv, sign);
0abe5247 1785 labelsv0daugh[ipr] = TMath::Abs(aodTrack->GetLabel());
ef87f0d6 1786 idv0daugh[ipr] = aodTrack->GetID();
0abe5247 1787 if (labelsv0daugh[ipr] == -1) isBkgV0 = kTRUE;
1788
1789 //Printf("v0 daughter %d has label %d", ipr, labelsv0daugh[ipr]);
ef87f0d6 1790
0abe5247 1791 AliKFParticle daughterKF(*aodTrack, pdg[ipr]); // we assume that the PDG is correct
1792 if (aodTrack->Charge() > 0) { // assigning positive and negative track to KF V0 for Armenteros-Podolanski plot
1793 positiveV0KF = daughterKF;
1794 }
1795 else {
1796 negativeV0KF = daughterKF;
1797 }
0abe5247 1798 }
ef87f0d6 1799
6be55058 1800 Double_t xn=0., xp=0.;//, dca;
ef87f0d6 1801 AliDebug(2, Form("bField = %f, esdv0Daugh1 = %p, esdv0Daugh2 = %p", fBField, esdv0Daugh1, esdv0Daugh2));
6be55058 1802 // dca = esdv0Daugh1->GetDCA(esdv0Daugh2, fBField, xn, xp);
ef87f0d6 1803
1804 AliExternalTrackParam tr1(*esdv0Daugh1);
1805 AliExternalTrackParam tr2(*esdv0Daugh2);
1806 tr1.PropagateTo(xn, fBField);
1807 tr2.PropagateTo(xp, fBField);
1808
1809 AliKFParticle daughterKF1(tr1, 211);
1810 AliKFParticle daughterKF2(tr2, 211);
1811 V0.AddDaughter(positiveV0KF);
1812 V0.AddDaughter(negativeV0KF);
1813 //V0.AddDaughter(daughterKF1);
1814 //V0.AddDaughter(daughterKF2);
1815
1816 delete esdv0Daugh1;
1817 delete esdv0Daugh2;
1818 esdv0Daugh1=0;
1819 esdv0Daugh2=0;
0abe5247 1820 // Checking the quality of the KF V0 vertex
1821 if( V0.GetNDF() < 1 ) {
1822 //Printf("Number of degrees of freedom < 1, continuing");
1823 return -1;
1824 }
ef87f0d6 1825 if( TMath::Sqrt(TMath::Abs(V0.GetChi2()/V0.GetNDF())) > fCutKFChi2NDF ) {
0abe5247 1826 //Printf("Chi2 per DOF too big, continuing");
1827 return -1;
1828 }
1829
1830 if(ftopoConstraint && nt > 0){
1831 for(Int_t ipr = 0; ipr < 2; ipr++){ // 0 is positive, 1 is negative
1832 AliAODTrack *aodTrack = (AliAODTrack*)v0part->GetDaughter(ipr);
1833 //* subtruct daughters from primary vertex
1834 if(!aodTrack->GetUsedForPrimVtxFit()) {
1835 //Printf("Track %d was not used for primary vertex, continuing", i);
1836 continue;
1837 }
1838 AliKFParticle daughterKF(*aodTrack, pdg[ipr]); // we assume that the PDG is correct
1839 primVtxCopy -= daughterKF;
1840 ntcheck--;
1841 }
1842 }
1843
0abe5247 1844 // Check V0 Chi^2 deviation from primary vertex // not needed for V0 for Lc decay!!
ef87f0d6 1845 /*
1846 if( V0.GetDeviationFromVertex( primVtxCopy ) < fCutKFDeviationFromVtxV0) {
0abe5247 1847 //Printf("Deviation from vertex too big, continuing");
0abe5247 1848 return -1;
1849 }
ef87f0d6 1850 */
0abe5247 1851
1852 //* Get V0 invariant mass
1853 Double_t massV0 = 999999, sigmaMassV0 = 999999;
1854 Int_t retMV0 = V0.GetMass( massV0, sigmaMassV0 );
1855 if( retMV0 ) {
1856 if (massV0 < 0) {
1857 codeKFV0 = 1; // Mass not ok
1858 if (sigmaMassV0 > 1e19) codeKFV0 = 5; // Mass and SigmaMass not ok
1859 }
1860 else if (sigmaMassV0 > 1e19) codeKFV0 = 2; // SigmaMass not ok
ef87f0d6 1861 }
0abe5247 1862 fHistoMassKFV0->Fill(massV0, sigmaMassV0);
1863
ef87f0d6 1864 if (massV0 < 0.4) Printf("\n\n>>>>>>>>>> Found the Funny V0 (mass = %f, sigma = %f, AOD mass = %f): labels of the tracks = %d, %d, id = %d and %d", massV0, sigmaMassV0, v0part->MassK0Short(), labelsv0daugh[0], labelsv0daugh[1], idv0daugh[0], idv0daugh[1]);
1865 if (massV0 > 0.55) Printf("\n\n>>>>>>>>>> Found the Funny V0 (mass = %f, , sigma = %f, AOD mass = %f): labels of the tracks = %d, %d, id = %d and %d", massV0, sigmaMassV0, v0part->MassK0Short(), labelsv0daugh[0], labelsv0daugh[1], idv0daugh[0], idv0daugh[1]);
1866
1867 Printf("Vertices: KF: x = %f, y = %f, z = %f", V0.GetX(), V0.GetY(), V0.GetZ());
1868 Printf("Vertices: AOD: x = %f, y = %f, z = %f", v0part->Xv(), v0part->Yv(), v0part->Zv());
1869
1870 //Printf("Got MC vtx for V0");
1871 if (fUseMCInfo && TMath::Abs(labelsv0daugh[0] - labelsv0daugh[1]) == 1) {
1872 AliAODMCParticle* tmpdaughv01 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelsv0daugh[0]));
1873 AliAODMCParticle* tmpdaughv02 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelsv0daugh[1]));
1874 if (!tmpdaughv01 && labelsv0daugh[0] > 0){
1875 AliDebug(2, "Could not access MC info for first daughter of V0, continuing");
1876 }
1877 if (!tmpdaughv02 && labelsv0daugh[1] > 0){
1878 AliDebug(2, "Could not access MC info for second daughter of V0, continuing");
1879 }
ff22c824 1880 if(tmpdaughv01){
1881 Double_t xPionMC = tmpdaughv01->Xv(); //Production vertex of Pion --> Where K0S decays
1882 Double_t yPionMC = tmpdaughv01->Yv();
1883 Double_t zPionMC = tmpdaughv01->Zv();
1884 //Printf("Got MC vtx for Pion");
1885 Printf("Vertices: MC: x = %f, y = %f, z = %f", xPionMC, yPionMC, zPionMC);
1886 }
ef87f0d6 1887 }
1888 else {
1889 Printf("Not a true V0");
1890 }
1891 //massV0=-1;//return -1;// !!!!
1892
1893 // now use what we just try with the bachelor, to build the Lc
1894
1895 // topological constraint
1896 nt = primVtxCopy.GetNContributors();
1897 ntcheck = nt;
1898
1899 Bool_t isMCokLc = kTRUE, isBkgLc = kFALSE;
1900 AliKFParticle Lc;
1901 Int_t labelsLcdaugh[2] = {-1, -1};
1902 labelsLcdaugh[0] = TMath::Abs(bach->GetLabel());
1903 labelsLcdaugh[1] = mcLabelV0;
1904
1905 if (bach->Charge() < 0) pdgLc[0] = -pdgLc[0];
1906 AliKFParticle daughterKFLc(*bach, pdgLc[0]);
1907 Lc.AddDaughter(daughterKFLc);
1908 TParticlePDG* particlePDG = TDatabasePDG::Instance()->GetParticle(310);
1909 Double_t massPDGK0S = particlePDG->Mass();
1910 V0.SetMassConstraint(massPDGK0S);
1911 Lc.AddDaughter(V0);
1912 if( Lc.GetNDF() < 1 ) {
1913 AliDebug(2, Form("Lc: Number of degrees of freedom < 1 (%d), continuing", Lc.GetNDF()));
1914 return -1;
1915 }
1916 if( TMath::Sqrt(TMath::Abs(Lc.GetChi2()/Lc.GetNDF())) > fCutKFChi2NDF) {
1917 AliDebug(2, Form("Lc: Chi2 per DOF too big, continuing (%f)", TMath::Sqrt(TMath::Abs(Lc.GetChi2()/Lc.GetNDF()))));
1918 return -1;
1919 }
1920
1921 if(ftopoConstraint && nt > 0){
1922 //* subtruct daughters from primary vertex
1923 if(!bach->GetUsedForPrimVtxFit()) {
1924 AliDebug(3, "Lc: Bachelor was not used for primary vertex, not subtracting it from primary vertex");
1925 }
1926 else{
1927 primVtxCopy -= daughterKFLc;
1928 ntcheck--;
1929 }
1930 /* the V0 was added above, so it is ok to remove it without checking
1931 if(!V0->GetUsedForPrimVtxFit()) {
1932 Printf("Lc: V0 was not used for primary vertex, continuing");
1933 continue;
1934 }
1935 */
1936 //primVtxCopy -= V0;
1937 //ntcheck--;
1938 }
1939
1940 // Check Lc Chi^2 deviation from primary vertex
1941 /*
1942 if( Lc.GetDeviationFromVertex( primVtxCopy ) > fCutKFDeviationFromVtx) {
1943 AliDebug(2, Form("Lc: Deviation from vertex too big, continuing (%f)", Lc.GetDeviationFromVertex( primVtxCopy )));
1944 return -1;
1945 }
1946
1947 if(ftopoConstraint){
1948 if(ntcheck>0) {
a3c87206 1949 // Add Lc to primary vertex to improve the primary vertex resolution
ef87f0d6 1950 primVtxCopy += Lc;
1951 Lc.SetProductionVertex(primVtxCopy);
1952 }
1953 }
1954 */
1955 //* Check chi^2
1956 if( TMath::Sqrt( TMath::Abs(Lc.GetChi2()/Lc.GetNDF())) > fCutKFChi2NDF) {
1957 AliDebug(2, Form("Lc: Final Chi2 per DOF too big, continuing (%f)", TMath::Sqrt( TMath::Abs(Lc.GetChi2()/Lc.GetNDF()))));
1958 return -1;
1959 }
1960
1961 if(ftopoConstraint){
1962 V0.SetProductionVertex(Lc);
1963 }
1964
1965 // After setting the vertex of the V0, getting/filling some info
1966
0abe5247 1967 //* Get V0 decayLength
1968 Double_t decayLengthV0 = 999999, sigmaDecayLengthV0 = 999999;
1969 Int_t retDLV0 = V0.GetDecayLength( decayLengthV0, sigmaDecayLengthV0 );
1970 if( retDLV0 ) {
1971 if (sigmaDecayLengthV0 > 1e19) {
1972 codeKFV0 = 3; // DecayLength not ok
1973 if (massV0 < 0) {
1974 codeKFV0 = 6; // DecayLength and Mass not ok
1975 if (sigmaMassV0 > 1e19) codeKFV0 = 11; // DecayLength and Mass and SigmaMass not ok
1976 }
1977 else if (sigmaMassV0 > 1e19) codeKFV0 = 8; // DecayLength and SigmaMass not ok
1978 }
1979 }
1980 fHistoDecayLengthKFV0->Fill(decayLengthV0, sigmaDecayLengthV0);
1981
1982 //* Get V0 life time
1983 Double_t lifeTimeV0 = 999999, sigmaLifeTimeV0 = 999999;
1984 Int_t retTLV0 = V0.GetLifeTime( lifeTimeV0, sigmaLifeTimeV0 );
1985 if( retTLV0 ) {
1986 if (sigmaLifeTimeV0 > 1e19) {
1987 codeKFV0 = 4; // LifeTime not ok
1988 if (sigmaDecayLengthV0 > 1e19) {
1989 codeKFV0 = 9; // LifeTime and DecayLength not ok
1990 if (massV0 < 0) {
1991 codeKFV0 = 14; // LifeTime and DecayLength and Mass not ok
1992 if (sigmaMassV0 > 1e19) codeKFV0 = 15; // LifeTime and DecayLength and Mass and SigmaMass not ok
1993 }
1994 else if (sigmaMassV0 > 1e19) codeKFV0 = 13; // LifeTime and DecayLength and SigmaMass not ok
1995 }
1996 else if (massV0 < 0) { // LifeTime and Mass and SigmaMass not ok
1997 codeKFV0 = 7; // LifeTime and Mass not ok
1998 if (sigmaMassV0 > 1e19) codeKFV0 = 12; // LifeTime and Mass and SigmaMass not ok
1999 }
2000 else if (sigmaMassV0 > 1e19) codeKFV0 = 10; // LifeTime and SigmaMass not ok
2001 }
2002 }
2003 fHistoLifeTimeKFV0->Fill(lifeTimeV0, sigmaLifeTimeV0);
2004
2005 if (codeKFV0 == -1) codeKFV0 = 0;
2006 fHistoKFV0->Fill(codeKFV0);
2007
2008 AliDebug(2, Form("V0: mass = %f, decay length = %f, life time = %f", massV0, decayLengthV0, lifeTimeV0 ));
2009
2010 fHistoMassV0All->Fill(massV0);
2011 fHistoDecayLengthV0All->Fill(decayLengthV0);
2012 fHistoLifeTimeV0All->Fill(lifeTimeV0);
2013
2014 Double_t qtAlphaV0[2] = {0., 0.};
2015 Double_t vtxV0KF[3] = {V0.GetX(), V0.GetY(), V0.GetZ()};
2016 positiveV0KF.TransportToPoint(vtxV0KF);
2017 negativeV0KF.TransportToPoint(vtxV0KF);
2018 V0.GetArmenterosPodolanski(positiveV0KF, negativeV0KF, qtAlphaV0);
2019 AliDebug(2, Form("Armenteros-Podolanski variables: alpha = %f, qt = %f", qtAlphaV0[1], qtAlphaV0[0]));
2020 fHistoArmenterosPodolanskiV0KF->Fill(qtAlphaV0[1], qtAlphaV0[0]);
2021 fHistoArmenterosPodolanskiV0AOD->Fill(v0part->AlphaV0(), v0part->PtArmV0());
2022 armPolKF[0] = qtAlphaV0[1];
2023 armPolKF[1] = qtAlphaV0[0];
2024
2025 // Checking MC info for V0
2026
2027 AliAODMCParticle *motherV0 = 0x0;
2028 AliAODMCParticle *daughv01 = 0x0;
2029 AliAODMCParticle *daughv02 = 0x0;
2030
2031 if (fUseMCInfo) {
2032 daughv01 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelsv0daugh[0]));
2033 daughv02 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelsv0daugh[1]));
2034 if (!daughv01 && labelsv0daugh[0] > 0){
2035 AliDebug(2, "Could not access MC info for first daughter of V0, continuing");
2036 isMCokV0 = kFALSE;
2037 }
2038 if (!daughv02 && labelsv0daugh[1] > 0){
2039 AliDebug(2, "Could not access MC info for second daughter of V0, continuing");
2040 isMCokV0 = kFALSE;
2041 }
2042 if (isMCokV0){
2043 if( daughv01->GetMother() == daughv02->GetMother() && daughv01->GetMother()>=0 ){
2044 AliDebug(3, Form("The mother has label %d", daughv01->GetMother()));
2045 motherV0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(daughv01->GetMother()));
2046 if( motherV0 && TMath::Abs(motherV0->GetPdgCode()) != 21 ){ // These are all V0 that are truly V0, not only K0S, but no gluons
2047 if( motherV0->GetNDaughters() == 2 ){
2048 fHistoMassV0True->Fill(massV0);
2049 fHistoDecayLengthV0True->Fill(decayLengthV0);
2050 fHistoLifeTimeV0True->Fill(lifeTimeV0);
2051 fHistoMassV0TrueFromAOD->Fill(v0part->MassK0Short());
2052 if (TMath::Abs(motherV0->GetPdgCode()) == 310){ // These are true V0 that are also K0S
2053 fHistoMassV0TrueK0S->Fill(massV0);
2054 fHistoDecayLengthV0TrueK0S->Fill(decayLengthV0);
2055 fHistoLifeTimeV0TrueK0S->Fill(lifeTimeV0);
2056 fHistoMassV0TrueK0SFromAOD->Fill(v0part->MassK0Short());
2057 }
2058 }
2059 AliDebug(2, Form("PDG V0 = %d, pi = %d, pj = %d, ndaughters = %d, mc mass = %f, reco mass = %f, v0 mass = %f", motherV0->GetPdgCode(), daughv01->GetPdgCode(), daughv02->GetPdgCode(), motherV0->GetNDaughters(), motherV0->GetCalcMass(), massV0, v0part->MassK0Short()));
2060 }
2061 else if (!motherV0){
2062 AliDebug(3, "could not access MC info for mother, continuing");
2063 }
2064 else {
2065 AliDebug(3, "MC mother is a gluon, continuing");
2066 }
2067 }
2068 else {
2069 AliDebug(3, "Background V0!");
2070 isBkgV0 = kTRUE;
2071 }
2072 }
2073 }
2074
2075 AliDebug(2, Form("isMCokV0 = %d, isBkgV0 = %d", (Int_t)isMCokV0, (Int_t)isBkgV0));
2076
ef87f0d6 2077 // Going back to Lc
0abe5247 2078
2079 //* Get Lc invariant mass
2080 Double_t massLc = 999999, sigmaMassLc= 999999;
2081 Int_t retMLc = Lc.GetMass( massLc, sigmaMassLc );
2082 if( retMLc ) {
2083 AliDebug(3, Form("----> Could not get mass (%e), and sigma(%e) for Lc, continuing", massLc, sigmaMassLc));
2084 if (massLc < 0) {
2085 codeKFLc = 1; // Mass not ok
2086 if (sigmaMassLc > 1e19) codeKFLc = 5; // Mass and SigmaMass not ok
2087 }
2088 else if (sigmaMassLc > 1e19) codeKFLc = 2; // SigmaMass not ok
2089 }
2090 fHistoMassKFLc->Fill(massLc, sigmaMassLc);
2091
2092 //* Get Lc Decay length
2093 Double_t decayLengthLc = 999999, sigmaDecayLengthLc = 999999;
2094 Int_t retDLLc = Lc.GetDecayLength( decayLengthLc, sigmaDecayLengthLc );
2095 if( retDLLc ) {
2096 AliDebug(3, "----> Lc: Could not get decay length, and sigma");
2097 if (sigmaDecayLengthLc > 1e19) {
2098 codeKFLc = 3; // DecayLength not ok
2099 if (massLc < 0) {
2100 codeKFLc = 6; // DecayLength and Mass not ok
2101 if (sigmaMassLc > 1e19) codeKFLc = 11; // DecayLength and Mass and SigmaMass not ok
2102 }
2103 else if (sigmaMassLc > 1e19) codeKFLc = 8; // DecayLength and SigmaMass not ok
2104 }
2105 }
2106 AliDebug(3, Form("retDLLc = %d, with decayLength = %f and error = %e", retDLLc, decayLengthLc, sigmaDecayLengthLc));
2107
2108 fHistoDecayLengthKFLc->Fill(decayLengthLc, sigmaDecayLengthLc);
2109
2110 //* Get Lc life time
2111 Double_t lifeTimeLc = 999999, sigmaLifeTimeLc = 999999;
2112 Int_t retTLLc = Lc.GetLifeTime( lifeTimeLc, sigmaLifeTimeLc );
2113 if( retTLLc ) {
2114 AliDebug(3, "----> Lc: Could not get lifeTime, and sigma");
2115 if (sigmaLifeTimeLc > 1e19) {
2116 codeKFLc = 4; // LifeTime not ok
2117 if (sigmaDecayLengthLc > 1e19) {
2118 codeKFLc = 9; // LifeTime and DecayLength not ok
2119 if (massLc < 0) {
2120 codeKFLc = 14; // LifeTime and DecayLength and Mass not ok
2121 if (sigmaMassLc > 1e19) codeKFLc = 15; // LifeTime and DecayLength and Mass and SigmaMass not ok
2122 }
2123 else if (sigmaMassLc > 1e19) codeKFLc = 13; // LifeTime and DecayLength and SigmaMass not ok
2124 }
2125 else if (massLc < 0) { // LifeTime and Mass and SigmaMass not ok
2126 codeKFLc = 7; // LifeTime and Mass not ok
2127 if (sigmaMassLc > 1e19) codeKFLc = 12; // LifeTime and Mass and SigmaMass not ok
2128 }
2129 else if (sigmaMassLc > 1e19) codeKFLc = 10; // LifeTime and SigmaMass not ok
2130 }
2131 }
2132
2133 fHistoLifeTimeKFLc->Fill(lifeTimeLc, sigmaLifeTimeLc);
2134
2135 AliDebug(2, Form("Lc: mass = %f (error = %e), decay length = %f (error = %e), life time = %f (error = %e) --> codeKFLc = %d", massLc, sigmaMassLc, decayLengthLc, sigmaDecayLengthLc, lifeTimeLc, sigmaLifeTimeLc, codeKFLc));
2136
2137 if (codeKFLc == -1) codeKFLc = 0;
2138 fHistoKFLc->Fill(codeKFLc);
2139
2140 fHistoKF->Fill(codeKFV0, codeKFLc);
2141
2142 // here we fill the histgrams for all the reconstructed KF vertices for the cascade
2143 fHistoMassLcAll->Fill(massLc);
2144 fHistoDecayLengthLcAll->Fill(decayLengthLc);
2145 fHistoLifeTimeLcAll->Fill(lifeTimeLc);
2146
2147 fHistoMassV0fromLcAll->Fill(massV0);
2148 fHistoDecayLengthV0fromLcAll->Fill(decayLengthV0);
2149 fHistoLifeTimeV0fromLcAll->Fill(lifeTimeV0);
2150
2151 Double_t xV0 = V0.GetX();
2152 Double_t yV0 = V0.GetY();
2153 Double_t zV0 = V0.GetZ();
2154
2155 Double_t xLc = Lc.GetX();
2156 Double_t yLc = Lc.GetY();
2157 Double_t zLc = Lc.GetZ();
2158
2159 Double_t xPrimVtx = primVtxCopy.GetX();
2160 Double_t yPrimVtx = primVtxCopy.GetY();
2161 Double_t zPrimVtx = primVtxCopy.GetZ();
2162
2163 Double_t distanceLcToPrimVtx = TMath::Sqrt((xPrimVtx - xLc) * (xPrimVtx - xLc) +
2164 (yPrimVtx - yLc) * (yPrimVtx - yLc) +
2165 (zPrimVtx - zLc) * (zPrimVtx - zLc));
2166
2167 Double_t distanceV0ToPrimVtx = TMath::Sqrt((xPrimVtx - xV0) * (xPrimVtx - xV0) +
2168 (yPrimVtx - yV0) * (yPrimVtx - yV0) +
2169 (zPrimVtx - zV0) * (zPrimVtx - zV0));
2170
2171 Double_t distanceV0ToLc = TMath::Sqrt((xLc - xV0)*(xLc - xV0) +
2172 (yLc - yV0)*(yLc - yV0) +
2173 (zLc - zV0)*(zLc - zV0));
2174
2175 //Printf("distanceLcToPrimVtx = %e, distanceV0ToPrimVtx= %f, distanceV0ToLc = %f", distanceLcToPrimVtx, distanceV0ToPrimVtx, distanceV0ToLc);
2176
2177 fHistoDistanceLcToPrimVtx->Fill(distanceLcToPrimVtx);
2178 fHistoDistanceV0ToPrimVtx->Fill(distanceV0ToPrimVtx);
2179 fHistoDistanceV0ToLc->Fill(distanceV0ToLc);
2180
2181 distances[0] = distanceLcToPrimVtx;
2182 distances[1] = distanceV0ToPrimVtx;
2183 distances[2] = distanceV0ToLc;
2184
2185 if (fUseMCInfo) {
2186
2187 AliAODMCParticle *daughv01Lc = 0x0;
2188 AliAODMCParticle *K0S = 0x0;
2189 AliAODMCParticle *daughv02Lc = 0x0;
2190
2191 if (labelsLcdaugh[0] >= 0) {
2192 // Printf("Getting Bachelor from label %d", labelsLcdaugh[1]);
2193 daughv01Lc = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelsLcdaugh[0]));
2194 if (!daughv01Lc){
2195 AliDebug(3, "Could not access MC info for first daughter of Lc");
2196 isMCokLc = kFALSE;
2197 }
2198 else {
2199 AliDebug(2, Form("The bachelor has label = %d", daughv01Lc->GetLabel()));
2200 if (TMath::Abs(daughv01Lc->GetPdgCode()) != 2212) isBkgLc = kTRUE;
2201 }
2202 }
2203 else { // The bachelor is a fake
2204 isBkgLc = kTRUE;
2205 }
2206
2207 if (labelsLcdaugh[1] >= 0) {
2208 //Printf("Getting K0S");
2209 K0S = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelsLcdaugh[1]));
2210 if (!K0S) {
2211 AliDebug(3, "Could not access MC info for second daughter of Lc");
2212 isMCokLc = kFALSE;
2213 }
2214 else{
2215 if (TMath::Abs(K0S->GetPdgCode()) != 310) isBkgLc = kTRUE;
2216 }
2217 }
2218 else{
2219 AliDebug(2, "The K0S is not true --> it does not have a label, continuing...");
2220 isBkgLc = kTRUE;
2221 }
2222
2223 if (!isBkgLc){ // so far, we only checked that the V0 and the bachelor are not fake, and in particular, we know that the V0 is a K0S since we used the MatchToMC method
2224 if (isMCokLc) { // We can then access its MC info, and it might then be that also the Lc is a true Lc
2225 Int_t iK0 = K0S->GetMother();
2226 if (iK0 < 0) {
2227 Printf("The K0S has no mother... IMPOSSIBLE"); // the K0S MUST have a mother!
2228 }
2229 else { // The K0S has a mother
2230 daughv02Lc = dynamic_cast<AliAODMCParticle*>(mcArray->At(iK0));
2231 if (!daughv02Lc){
2232 AliDebug(3, "Could not access MC info for second daughter of Lc");
2233 }
2234 else { // we can access safely the K0S mother in the MC
ae48be26 2235 if( daughv01Lc && (daughv01Lc->GetMother() == daughv02Lc->GetMother()) && (daughv01Lc->GetMother()>=0) ){ // This is a true cascade! bachelor and V0 come from the same mother
0abe5247 2236 //Printf("Lc: The mother has label %d", daughv01Lc->GetMother());
2237 AliAODMCParticle *motherLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(daughv01Lc->GetMother()));
2238 Int_t pdgMum = 0, pdgBach = 0, pdgV0 = 0;
2239 if (motherLc) pdgMum = motherLc->GetPdgCode();
2240 if (daughv01Lc) pdgBach = daughv01Lc->GetPdgCode();
2241 if (daughv02Lc) pdgV0 = daughv02Lc->GetPdgCode();
2242 AliDebug(2, Form("pdgMum = %d, pdgBach = %d, pdgV0 = %d", pdgMum, pdgBach, pdgV0));
2243
2244 if( motherLc && TMath::Abs(motherLc->GetPdgCode()) != 21 ){ // These are true cascades, we don't know yet if they are Lc
2245 fHistoMassLcTrue->Fill(massLc);
2246 fHistoDecayLengthLcTrue->Fill(decayLengthLc);
2247 fHistoLifeTimeLcTrue->Fill(lifeTimeLc);
2248 fHistoMassLcTrueFromAOD->Fill(cascade->InvMassLctoK0sP());
2249
2250 fHistoMassV0fromLcTrue->Fill(massV0);
2251 fHistoDecayLengthV0fromLcTrue->Fill(decayLengthV0);
2252 fHistoLifeTimeV0fromLcTrue->Fill(lifeTimeV0);
2253
2254 if (TMath::Abs(motherLc->GetPdgCode()) == 4122 && TMath::Abs(motherV0->GetPdgCode()) == 310 && TMath::Abs(daughv01Lc->GetPdgCode()) == 2212){ // This is Lc --> K0S + p (the check on the PDG code of the V0 is useless, since we used MathcToMC with it, but fine...
2255 AliDebug(2, Form("IT IS SIGNAL!!! with label = %d", motherLc->GetLabel()));
2256
2257 fHistoArmenterosPodolanskiV0KFSgn->Fill(qtAlphaV0[1], qtAlphaV0[0]);
2258 fHistoArmenterosPodolanskiV0AODSgn->Fill(v0part->AlphaV0(), v0part->PtArmV0());
2259
2260 fHistoDistanceLcToPrimVtxSgn->Fill(distanceLcToPrimVtx);
2261 fHistoDistanceV0ToPrimVtxSgn->Fill(distanceV0ToPrimVtx);
2262 fHistoDistanceV0ToLcSgn->Fill(distanceV0ToLc);
2263
2264 fHistoMassLcSgn->Fill(massLc);
2265 fHistoMassLcSgnFromAOD->Fill(cascade->InvMassLctoK0sP());
2266
2267 fHistoDecayLengthLcSgn->Fill(decayLengthLc);
2268 fHistoLifeTimeLcSgn->Fill(lifeTimeLc);
2269
2270 fHistoMassV0fromLcSgn->Fill(massV0);
2271 fHistoDecayLengthV0fromLcSgn->Fill(decayLengthV0);
2272 fHistoLifeTimeV0fromLcSgn->Fill(lifeTimeV0);
2273 }
2274 else {
2275 isBkgLc = kTRUE; // it is not a Lc, since the pdg != 4122
2276 }
2277
2278 // if we got here, we can compare with MC information; this part is done also for cases where the cascade is not a Lc!
2279 // First, we compare the vtx of the Lc
2280 Double_t xLcMC = motherLc->Xv();
2281 Double_t yLcMC = motherLc->Yv();
2282 Double_t zLcMC = motherLc->Zv();
2283 //Printf("Got MC vtx for Lc");
2284 Double_t xProtonMC = daughv01Lc->Xv();
2285 Double_t yProtonMC = daughv01Lc->Yv();
2286 Double_t zProtonMC = daughv01Lc->Zv();
2287 //Printf("Got MC vtx for Proton");
2288
2289 Double_t vtxLcResidualToPrimVtx = TMath::Sqrt((xLcMC - xProtonMC) * (xLcMC - xProtonMC) +
2290 (yLcMC - yProtonMC) * (yLcMC - yProtonMC) +
2291 (zLcMC - zProtonMC) * (zLcMC - zProtonMC)) - distanceLcToPrimVtx;
2292
2293 // Then, we compare the vtx of the K0s
2294 Double_t xV0MC = motherV0->Xv(); //Production vertex of K0S --> Where Lc decays
2295 Double_t yV0MC = motherV0->Yv();
2296 Double_t zV0MC = motherV0->Zv();
ef87f0d6 2297
0abe5247 2298 //Printf("Got MC vtx for V0");
2299 Double_t xPionMC = daughv01->Xv(); //Production vertex of Pion --> Where K0S decays
2300 Double_t yPionMC = daughv01->Yv();
2301 Double_t zPionMC = daughv01->Zv();
2302 //Printf("Got MC vtx for Pion");
ef87f0d6 2303 Printf("Vertices: MC: x = %f, y = %f, z = %f", xPionMC, yPionMC, zPionMC);
0abe5247 2304
2305 Double_t vtxV0ResidualToLc = TMath::Sqrt((xV0MC - xPionMC) * (xV0MC - xPionMC) +
2306 (yV0MC - yPionMC) * (yV0MC - yPionMC) +
2307 (zV0MC - zPionMC) * (zV0MC - zPionMC)) - distanceV0ToLc;
2308
2309 // Then, the K0S vertex wrt primary vertex
2310
2311 Double_t vtxV0ResidualToPrimVtx = TMath::Sqrt((xPionMC - xLcMC) * (xPionMC - xLcMC) +
2312 (yPionMC - yLcMC) * (yPionMC - yLcMC) +
2313 (zPionMC - zLcMC) * (zPionMC - zLcMC)) - distanceV0ToPrimVtx;
2314
2315 fHistoVtxLcResidualToPrimVtx->Fill(vtxLcResidualToPrimVtx);
2316 fHistoVtxV0ResidualToLc->Fill(vtxV0ResidualToLc);
2317 fHistoVtxV0ResidualToPrimVtx->Fill(vtxV0ResidualToPrimVtx);
2318
2319 } //closing: if( motherLc && TMath::Abs(motherLc->GetPdgCode()) != 21 )
2320 else if (!motherLc){
2321 AliDebug(2, "We could not access MC info for Lc mother, so we did nothing");
2322 }
2323 else {
2324 AliDebug(2, "MC Lc mother is a gluon, so we did nothing");
2325 }
2326 } //closing: if( daughv01Lc->GetMother() == daughv02Lc->GetMother() && daughv01Lc->GetMother()>=0 )
2327 else {
2328 isBkgLc = kTRUE; // it cannot be a Lc, since the daughters do not have the same mother
2329 }
2330 } // closing: else { // we can access safely the K0S mother in the MC
2331 } // closing: else { // The K0S has a mother
2332 } // closing isMCLcok
2333 } // closing !isBkgLc
2334 } // closing fUseMCInfo
2335
2336 //Printf("retMV0 = %d, retMLc = %d", retMV0, retMLc);
2337 if ( retMV0 == 0 && retMLc == 0){
2338 V0KF[0] = massV0;
2339 errV0KF[0] = sigmaMassV0;
2340 V0KF[1] = decayLengthV0;
2341 errV0KF[1] = sigmaDecayLengthV0;
2342 V0KF[2] = lifeTimeV0;
2343 errV0KF[2] = sigmaLifeTimeV0;
2344 LcKF[0] = massLc;
2345 errLcKF[0] = sigmaMassLc;
2346 LcKF[1] = decayLengthLc;
2347 errLcKF[1] = sigmaDecayLengthLc;
2348 LcKF[2] = lifeTimeLc;
2349 errLcKF[2] = sigmaLifeTimeLc;
2350 }
2351
2352 return 1;
2353
2354}
2355//________________________________________________________________________
2356AliAnalysisTaskSELc2V0bachelorTMVA::EBachelor AliAnalysisTaskSELc2V0bachelorTMVA::CheckBachelor( AliAODRecoCascadeHF *part,
2357 AliAODTrack* bachelor,
2358 TClonesArray *mcArray ){
2359
2360 //Printf("In CheckBachelor");
2361
2362 // function to check if the bachelor is a p, if it is a p but not from Lc
2363 // to be filled for background candidates
2364
2365 Int_t label = bachelor->GetLabel();
2366 if (label == -1) {
2367 return kBachFake;
2368 }
2369
2370 AliAODMCParticle *mcpart = dynamic_cast<AliAODMCParticle*>(mcArray->At(TMath::Abs(label)));
2371 if (!mcpart) return kBachInvalid;
2372 Int_t pdg = mcpart->PdgCode();
2373 if (TMath::Abs(pdg) != 2212) {
2374 AliDebug(2, Form("Bachelor is not a p, but a particle with pdg code = %d", pdg));
2375 return kBachNoProton;
2376 }
2377 else { // it is a proton
2378 //Int_t labelLc = part->GetLabel();
2379 Int_t labelLc = FindLcLabel(part, mcArray);
2380 //Printf(">>>>>>>>>>>>> label for Lc = %d", labelLc);
2381 Int_t bachelorMotherLabel = mcpart->GetMother();
2382 //Printf(">>>>>>>>>>>>> label for bachelorMother = %d", bachelorMotherLabel);
2383 if (bachelorMotherLabel == -1) {
2384 return kBachPrimary;
2385 }
2386 AliAODMCParticle *bachelorMother = dynamic_cast<AliAODMCParticle*>(mcArray->At(bachelorMotherLabel));
2387 if (!bachelorMother) return kBachInvalid;
2388 Int_t pdgMother = bachelorMother->PdgCode();
2389 if (TMath::Abs(pdgMother) != 4122) {
2390 AliDebug(2, Form("The proton does not come from a Lc, but from a particle with pdgcode = %d", pdgMother));
2391 return kBachNoLambdaMother;
2392 }
2393 else { // it comes from Lc
2394 if (labelLc != bachelorMotherLabel){
2395 //AliInfo(Form("The proton comes from a Lc, but it is not the candidate we are analyzing (label Lc = %d, label p mother = %d", labelLc, bachelorMotherLabel));
2396 AliDebug(2, Form("The proton comes from a Lc, but it is not the candidate we are analyzing"));
2397 return kBachDifferentLambdaMother;
2398 }
2399 else { // it comes from the correct Lc
2400 AliDebug(2, Form("The proton comes from a Lc, and it is EXACTLY the candidate we are analyzing"));
2401 return kBachCorrectLambdaMother;
2402 }
2403 }
2404 }
2405
2406 return kBachInvalid;
2407
2408}
2409
2410//________________________________________________________________________
2411AliAnalysisTaskSELc2V0bachelorTMVA::EK0S AliAnalysisTaskSELc2V0bachelorTMVA::CheckK0S( AliAODRecoCascadeHF *part,
2412 AliAODv0* v0part,
2413 //AliAODTrack* v0part,
2414 TClonesArray *mcArray ){
2415
2416 // function to check if the K0Spart is a p, if it is a p but not from Lc
2417 // to be filled for background candidates
2418
2419 //Printf(" CheckK0S");
2420
2421 //Int_t labelMatchToMC = v0part->MatchToMC(310, mcArray);
2422 //Int_t label = v0part->GetLabel();
2423 Int_t labelFind = FindV0Label(v0part, mcArray);
2424 //Printf("\n >>>>>>>>>>>>> label for V0 = %d, from MatchToMC = %d, from Find = %d", label, labelMatchToMC, labelFind);
2425 if (labelFind == -1) {
2426 return kK0SFake;
2427 }
2428
2429 AliAODMCParticle *mcpart = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelFind));
2430 if (!mcpart) return kK0SInvalid;
2431 Int_t pdg = mcpart->PdgCode();
2432 if (TMath::Abs(pdg) != 310) {
2433 AliDebug(2, Form("K0Spart is not a K0S, but a particle with pdg code = %d", pdg));
2434 //AliInfo(Form("K0Spart is not a K0S, but a particle with pdg code = %d", pdg));
2435 return kK0SNoK0S;
2436 }
2437 else { // it is a K0S
2438 //Int_t labelLc = part->GetLabel();
2439 Int_t labelLc = FindLcLabel(part, mcArray);
2440 Int_t K0SpartMotherLabel = mcpart->GetMother();
2441 if (K0SpartMotherLabel == -1) {
2442 return kK0SWithoutMother;
2443 }
2444 AliAODMCParticle *K0SpartMother = dynamic_cast<AliAODMCParticle*>(mcArray->At(K0SpartMotherLabel)); // should be a K0 in case of signal
2445 if (!K0SpartMother) return kK0SInvalid;
2446 Int_t pdgMotherK0S = K0SpartMother->PdgCode();
2447 if (TMath::Abs(pdgMotherK0S) != 311) {
2448 AliDebug(2, Form("The K0S does not come from a K0, but from a particle with pdgcode = %d", pdgMotherK0S));
2449 //AliInfo(Form("The K0S does not come from a K0, but from a particle with pdgcode = %d", pdgMotherK0S));
2450 return kK0SNotFromK0;
2451 }
2452 else { // the K0S comes from a K0
2453 Int_t K0MotherLabel = K0SpartMother->GetMother(); // mother of K0 --> Lc in case of signal
2454 if (K0MotherLabel == -1) {
2455 return kK0Primary;
2456 }
2457 AliAODMCParticle *K0Mother = dynamic_cast<AliAODMCParticle*>(mcArray->At(K0MotherLabel));
2458 if (!K0Mother) return kK0SInvalid;
2459 Int_t pdgK0Mother = K0Mother->PdgCode();
2460 if (TMath::Abs(pdgK0Mother) != 4122) { // the K0 does not come from a Lc
2461 AliDebug(2, Form("The K0 does not come from a Lc, but from a particle with pdgcode = %d", pdgK0Mother));
2462 //AliInfo(Form("The K0 does not come from a Lc, but from a particle with pdgcode = %d", pdgK0Mother));
2463 return kK0NoLambdaMother;
2464 }
2465 else { // the K0 comes from Lc
2466 if (labelLc != K0MotherLabel){ // The K0 comes from a different Lc
2467 AliDebug(2, Form("The K0S comes from a Lc, but it is not the candidate we are analyzing"));
2468 //AliInfo(Form("The K0S comes from a Lc, but it is not the candidate we are analyzing"));
2469 return kK0DifferentLambdaMother;
2470 }
2471 else { // it comes from the correct Lc
2472 AliDebug(2, Form("The K0S comes from a Lc, and it is EXACTLY the candidate we are analyzing"));
2473 //AliInfo(Form("The K0S comes from a Lc, and it is EXACTLY the candidate we are analyzing"));
2474 return kK0CorrectLambdaMother;
2475 }
2476 }
2477 }
2478 }
2479
2480 return kK0SInvalid;
2481
2482}
2483
2484//----------------------------------------------------------------------------
2485Int_t AliAnalysisTaskSELc2V0bachelorTMVA::FindV0Label(AliAODRecoDecay* v0part, TClonesArray *mcArray) const
2486{
2487
2488 //Printf(" FindV0Label");
2489
2490 // finding the label of teh V0; inspired from AliAODRecoDecay::MatchToMC
2491
2492 Int_t labMother[2]={-1, -1};
2493 AliAODMCParticle *part=0;
2494 AliAODMCParticle *mother=0;
2495 Int_t dgLabels = -1;
2496
2497 Int_t ndg = v0part->GetNDaughters();
2498 if (ndg != 2) {
2499 AliFatal(Form("IMPOSSIBLE!! There are %d daughters, but they should be 2!", ndg));
2500 }
2501
2502 for(Int_t i = 0; i < 2; i++) {
2503 AliAODTrack *trk = (AliAODTrack*)v0part->GetDaughter(i);
2504 dgLabels = trk->GetLabel();
2505 if (dgLabels == -1) {
2506 //printf("daughter with negative label %d\n", dgLabels);
2507 return -1;
2508 }
2509 part = (AliAODMCParticle*)mcArray->At(TMath::Abs(dgLabels));
2510 if (!part) {
2511 //printf("no MC particle\n");
2512 return -1;
2513 }
2514 labMother[i] = part->GetMother();
2515 if (labMother[i] != -1){
2516 mother = (AliAODMCParticle*)mcArray->At(labMother[i]);
2517 if(!mother) {
2518 //printf("no MC mother particle\n");
2519 return -1;
2520 }
2521 }
2522 else {
2523 return -1;
2524 }
2525 }
2526
2527 if (labMother[0] == labMother[1]) return labMother[0];
2528
2529 else return -1;
2530
2531}
2532
2533//----------------------------------------------------------------------------
2534Int_t AliAnalysisTaskSELc2V0bachelorTMVA::FindLcLabel(AliAODRecoCascadeHF* cascade, TClonesArray *mcArray) const
2535{
2536
2537 // finding the label of teh V0; inspired from AliAODRecoDecay::MatchToMC
2538
2539 //Printf(" FindLcLabel");
2540
2541 AliAODMCParticle *part=0;
2542 AliAODMCParticle *mother=0;
2543 AliAODMCParticle *grandmother=0;
2544 Int_t dgLabels = -1;
2545
2546 Int_t ndg = cascade->GetNDaughters();
2547 if (ndg != 2) {
2548 AliFatal(Form("IMPOSSIBLE!! There are %d daughters, but they should be 2!", ndg));
2549 }
2550
2551 // daughter 0 --> bachelor, daughter 1 --> V0
2552 AliAODTrack *trk = (AliAODTrack*)cascade->GetDaughter(0); // bachelor
2553 dgLabels = trk->GetLabel();
2554 if (dgLabels == -1 ) {
2555 //printf("daughter with negative label %d\n", dgLabels);
2556 return -1;
2557 }
2558 part = (AliAODMCParticle*)mcArray->At(TMath::Abs(dgLabels));
2559 if (!part) {
2560 //printf("no MC particle\n");
2561 return -1;
2562 }
2563 Int_t labMotherBach = part->GetMother();
2564 if (labMotherBach == -1){
2565 return -1;
2566 }
2567 mother = (AliAODMCParticle*)mcArray->At(labMotherBach);
2568 if(!mother) {
2569 //printf("no MC mother particle\n");
2570 return -1;
2571 }
2572
2573 AliAODv0 *v0 = (AliAODv0*)cascade->GetDaughter(1); // V0
2574 dgLabels = FindV0Label(v0, mcArray);
2575 if (dgLabels == -1 ) {
2576 //printf("daughter with negative label (v0 was a fake) %d\n", dgLabels);
2577 return -1;
2578 }
2579 part = (AliAODMCParticle*)mcArray->At(TMath::Abs(dgLabels));
2580 if (!part) {
2581 //printf("no MC particle\n");
2582 return -1;
2583 }
2584 Int_t labMotherv0 = part->GetMother();
2585 if (labMotherv0 == -1){
2586 return -1;
2587 }
2588 mother = (AliAODMCParticle*)mcArray->At(labMotherv0);
2589 if(!mother) {
2590 //printf("no MC mother particle\n");
2591 return -1;
2592 }
2593 Int_t labGrandMotherv0 = mother->GetMother();
2594 if (labGrandMotherv0 == -1){
2595 return -1;
2596 }
2597 grandmother = (AliAODMCParticle*)mcArray->At(labGrandMotherv0);
2598 if(!grandmother) {
2599 //printf("no MC mother particle\n");
2600 return -1;
2601 }
2602
2603 //Printf("labMotherBach = %d, labMotherv0 = %d, labGrandMotherv0 = %d", labMotherBach, labMotherv0, labGrandMotherv0);
2604 if (labGrandMotherv0 == labMotherBach) return labMotherBach;
2605
2606 else return -1;
2607
2608}
0abe5247 2609
2610
2611
2612
2613
2614