AliPhysics  b11e70a (b11e70a)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisVertexingHF.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 /* $Id$ */
17 
18 //----------------------------------------------------------------------------
19 // Implementation of the heavy-flavour vertexing analysis class
20 // Candidates are stored in the AOD as objects deriving from AliAODRecoDecay.
21 // To be used as a task of AliAnalysisManager by means of the interface
22 // class AliAnalysisTaskSEVertexingHF.
23 // An example of usage in the macro AliAnalysisTaskSEVertexingHFTest.C.
24 //
25 // Contact: andrea.dainese@pd.infn.it
26 // Contributors: E.Bruna, G.E.Bruno, A.Dainese, C.Di Gliglio,
27 // F.Prino, R.Romita, X.M.Zhang
28 //----------------------------------------------------------------------------
29 #include <Riostream.h>
30 #include <TFile.h>
31 #include <TDatabasePDG.h>
32 #include <TString.h>
33 #include <TList.h>
34 #include <TProcessID.h>
35 #include "AliLog.h"
36 #include "AliVEvent.h"
37 #include "AliVVertex.h"
38 #include "AliVTrack.h"
39 #include "AliVertexerTracks.h"
40 #include "AliKFVertex.h"
41 #include "AliESDEvent.h"
42 #include "AliESDVertex.h"
43 #include "AliExternalTrackParam.h"
44 #include "AliNeutralTrackParam.h"
45 #include "AliESDtrack.h"
46 #include "AliESDtrackCuts.h"
47 #include "AliAODEvent.h"
48 #include "AliPIDResponse.h"
49 #include "AliAODRecoDecay.h"
50 #include "AliAODRecoDecayHF.h"
54 #include "AliAODRecoCascadeHF.h"
55 #include "AliRDHFCutsD0toKpi.h"
56 #include "AliRDHFCutsJpsitoee.h"
58 #include "AliRDHFCutsDstoKKpi.h"
59 #include "AliRDHFCutsLctopKpi.h"
60 #include "AliRDHFCutsLctoV0.h"
61 #include "AliRDHFCutsD0toKpipipi.h"
63 #include "AliAnalysisFilter.h"
64 #include "AliAnalysisVertexingHF.h"
65 #include "AliMixedEvent.h"
66 #include "AliESDv0.h"
67 #include "AliAODv0.h"
68 #include "AliCodeTimer.h"
69 #include "AliMultSelection.h"
70 #include <cstring>
71 
75 
76 //----------------------------------------------------------------------------
78 fInputAOD(kFALSE),
79 fAODMapSize(0),
80 fAODMap(0),
81 fVertexerTracks(0x0),
82 fBzkG(0.),
83 fSecVtxWithKF(kFALSE),
84 fRecoPrimVtxSkippingTrks(kFALSE),
85 fRmTrksFromPrimVtx(kFALSE),
86 fV1(0x0),
87 fD0toKpi(kTRUE),
88 fJPSItoEle(kTRUE),
89 f3Prong(kTRUE),
90 f4Prong(kTRUE),
91 fDstar(kTRUE),
92 fCascades(kTRUE),
93 fLikeSign(kFALSE),
94 fLikeSign3prong(kFALSE),
95 fMixEvent(kFALSE),
96 fPidResponse(0x0),
97 fUseKaonPIDfor3Prong(kFALSE),
98 fUsePIDforLc(0),
99 fUsePIDforLc2V0(0),
100 fUseKaonPIDforDs(kFALSE),
101 fUseTPCPID(kFALSE),
102 fUseTOFPID(kFALSE),
103 fUseTPCPIDOnlyIfNoTOF(kFALSE),
104 fMaxMomForTPCPid(1.),
105 fUsePidTag(kFALSE),
106 fnSigmaTPCPionLow(5.),
107 fnSigmaTPCPionHi(5.),
108 fnSigmaTOFPionLow(5.),
109 fnSigmaTOFPionHi(5.),
110 fnSigmaTPCKaonLow(5.),
111 fnSigmaTPCKaonHi(5.),
112 fnSigmaTOFKaonLow(5.),
113 fnSigmaTOFKaonHi(5.),
114 fnSigmaTPCProtonLow(5.),
115 fnSigmaTPCProtonHi(5.),
116 fnSigmaTOFProtonLow(5.),
117 fnSigmaTOFProtonHi(5.),
118 fMaxCentPercentileForTightCuts(-9999),
119 fTrackFilter(0x0),
120 fTrackFilter2prongCentral(0x0),
121 fTrackFilter3prongCentral(0x0),
122 fTrackFilterSoftPi(0x0),
123 fCutsD0toKpi(0x0),
124 fCutsJpsitoee(0x0),
125 fCutsDplustoKpipi(0x0),
126 fCutsDstoKKpi(0x0),
127 fCutsLctopKpi(0x0),
128 fCutsLctoV0(0x0),
129 fCutsD0toKpipipi(0x0),
130 fCutsDStartoKpipi(0x0),
131 fListOfCuts(0x0),
132 fFindVertexForDstar(kTRUE),
133 fFindVertexForCascades(kTRUE),
134 fV0TypeForCascadeVertex(0),
135 fMassCutBeforeVertexing(kFALSE),
136 fMassCalc2(0),
137 fMassCalc3(0),
138 fMassCalc4(0),
139 fOKInvMassD0(kFALSE),
140 fOKInvMassJpsi(kFALSE),
141 fOKInvMassDplus(kFALSE),
142 fOKInvMassDs(kFALSE),
143 fOKInvMassLc(kFALSE),
144 fOKInvMassDstar(kFALSE),
145 fOKInvMassD0to4p(kFALSE),
146 fOKInvMassLctoV0(kFALSE),
147 fnTrksTotal(0),
148 fnSeleTrksTotal(0),
149 fMakeReducedRHF(kFALSE),
150 fMassDzero(0.),
151 fMassDplus(0.),
152 fMassDs(0.),
153 fMassLambdaC(0.),
154 fMassDstar(0.),
155 fMassJpsi(0.)
156 {
158 
159  Double_t d02[2]={0.,0.};
160  Double_t d03[3]={0.,0.,0.};
161  Double_t d04[4]={0.,0.,0.,0.};
162  fMassCalc2 = new AliAODRecoDecay(0x0,2,0,d02);
163  fMassCalc3 = new AliAODRecoDecay(0x0,3,1,d03);
164  fMassCalc4 = new AliAODRecoDecay(0x0,4,0,d04);
165  SetMasses();
166 }
167 //--------------------------------------------------------------------------
169 TNamed(source),
170 fInputAOD(source.fInputAOD),
171 fAODMapSize(source.fAODMapSize),
172 fAODMap(source.fAODMap),
173 fVertexerTracks(source.fVertexerTracks),
174 fBzkG(source.fBzkG),
175 fSecVtxWithKF(source.fSecVtxWithKF),
176 fRecoPrimVtxSkippingTrks(source.fRecoPrimVtxSkippingTrks),
177 fRmTrksFromPrimVtx(source.fRmTrksFromPrimVtx),
178 fV1(source.fV1),
179 fD0toKpi(source.fD0toKpi),
180 fJPSItoEle(source.fJPSItoEle),
181 f3Prong(source.f3Prong),
182 f4Prong(source.f4Prong),
183 fDstar(source.fDstar),
184 fCascades(source.fCascades),
185 fLikeSign(source.fLikeSign),
186 fLikeSign3prong(source.fLikeSign3prong),
187 fMixEvent(source.fMixEvent),
188 fPidResponse(source.fPidResponse),
189 fUseKaonPIDfor3Prong(source.fUseKaonPIDfor3Prong),
190 fUsePIDforLc(source.fUsePIDforLc),
191 fUsePIDforLc2V0(source.fUsePIDforLc2V0),
192 fUseKaonPIDforDs(source.fUseKaonPIDforDs),
193 fUseTPCPID(source.fUseTPCPID),
194 fUseTOFPID(source.fUseTOFPID),
195 fUseTPCPIDOnlyIfNoTOF(source.fUseTPCPIDOnlyIfNoTOF),
196 fMaxMomForTPCPid(source.fMaxMomForTPCPid),
197 fUsePidTag(source.fUsePidTag),
198 fnSigmaTPCPionLow(source.fnSigmaTPCPionLow),
199 fnSigmaTPCPionHi(source.fnSigmaTPCPionHi),
200 fnSigmaTOFPionLow(source.fnSigmaTOFPionLow),
201 fnSigmaTOFPionHi(source.fnSigmaTOFPionHi),
202 fnSigmaTPCKaonLow(source.fnSigmaTPCKaonLow),
203 fnSigmaTPCKaonHi(source.fnSigmaTPCKaonHi),
204 fnSigmaTOFKaonLow(source.fnSigmaTOFKaonLow),
205 fnSigmaTOFKaonHi(source.fnSigmaTOFKaonHi),
206 fnSigmaTPCProtonLow(source.fnSigmaTPCProtonLow),
207 fnSigmaTPCProtonHi(source.fnSigmaTPCProtonHi),
208 fnSigmaTOFProtonLow(source.fnSigmaTOFProtonLow),
209 fnSigmaTOFProtonHi(source.fnSigmaTOFProtonHi),
210 fMaxCentPercentileForTightCuts(source.fMaxCentPercentileForTightCuts),
211 fTrackFilter(source.fTrackFilter),
212 fTrackFilter2prongCentral(source.fTrackFilter2prongCentral),
213 fTrackFilter3prongCentral(source.fTrackFilter3prongCentral),
214 fTrackFilterSoftPi(source.fTrackFilterSoftPi),
215 fCutsD0toKpi(source.fCutsD0toKpi),
216 fCutsJpsitoee(source.fCutsJpsitoee),
217 fCutsDplustoKpipi(source.fCutsDplustoKpipi),
218 fCutsDstoKKpi(source.fCutsDstoKKpi),
219 fCutsLctopKpi(source.fCutsLctopKpi),
220 fCutsLctoV0(source.fCutsLctoV0),
221 fCutsD0toKpipipi(source.fCutsD0toKpipipi),
222 fCutsDStartoKpipi(source.fCutsDStartoKpipi),
223 fListOfCuts(source.fListOfCuts),
224 fFindVertexForDstar(source.fFindVertexForDstar),
225 fFindVertexForCascades(source.fFindVertexForCascades),
226 fV0TypeForCascadeVertex(source.fV0TypeForCascadeVertex),
227 fMassCutBeforeVertexing(source.fMassCutBeforeVertexing),
228 fMassCalc2(source.fMassCalc2),
229 fMassCalc3(source.fMassCalc3),
230 fMassCalc4(source.fMassCalc4),
231 fOKInvMassD0(source.fOKInvMassD0),
232 fOKInvMassJpsi(source.fOKInvMassJpsi),
233 fOKInvMassDplus(source.fOKInvMassDplus),
234 fOKInvMassDs(source.fOKInvMassDs),
235 fOKInvMassLc(source.fOKInvMassLc),
236 fOKInvMassDstar(source.fOKInvMassDstar),
237 fOKInvMassD0to4p(source.fOKInvMassD0to4p),
238 fOKInvMassLctoV0(source.fOKInvMassLctoV0),
239 fnTrksTotal(0),
240 fnSeleTrksTotal(0),
241 fMakeReducedRHF(kFALSE),
242 fMassDzero(source.fMassDzero),
243 fMassDplus(source.fMassDplus),
244 fMassDs(source.fMassDs),
245 fMassLambdaC(source.fMassLambdaC),
246 fMassDstar(source.fMassDstar),
247 fMassJpsi(source.fMassJpsi)
248 {
252 }
253 //--------------------------------------------------------------------------
255 {
259  if(&source == this) return *this;
260  fInputAOD = source.fInputAOD;
261  fAODMapSize = source.fAODMapSize;
263  fBzkG = source.fBzkG;
264  fSecVtxWithKF = source.fSecVtxWithKF;
267  fV1 = source.fV1;
268  fD0toKpi = source.fD0toKpi;
269  fJPSItoEle = source.fJPSItoEle;
270  f3Prong = source.f3Prong;
271  f4Prong = source.f4Prong;
272  fDstar = source.fDstar;
273  fCascades = source.fCascades;
274  fLikeSign = source.fLikeSign;
276  fMixEvent = source.fMixEvent;
277  fPidResponse = source.fPidResponse;
279  fUsePIDforLc = source.fUsePIDforLc;
282  fUseTPCPID = source.fUseTPCPID;
283  fUseTOFPID = source.fUseTOFPID;
286  fUsePidTag = source.fUsePidTag;
302  fTrackFilter = source.fTrackFilter;
306  fCutsD0toKpi = source.fCutsD0toKpi;
307  fCutsJpsitoee = source.fCutsJpsitoee;
309  fCutsDstoKKpi = source.fCutsDstoKKpi;
310  fCutsLctopKpi = source.fCutsLctopKpi;
311  fCutsLctoV0 = source.fCutsLctoV0;
314  fListOfCuts = source.fListOfCuts;
319  fMassCalc2 = source.fMassCalc2;
320  fMassCalc3 = source.fMassCalc3;
321  fMassCalc4 = source.fMassCalc4;
322  fOKInvMassD0 = source.fOKInvMassD0;
325  fOKInvMassDs = source.fOKInvMassDs;
326  fOKInvMassLc = source.fOKInvMassLc;
330  fMassDzero = source.fMassDzero;
331  fMassDplus = source.fMassDplus;
332  fMassDs = source.fMassDs;
333  fMassLambdaC = source.fMassLambdaC;
334  fMassDstar = source.fMassDstar;
335  fMassJpsi = source.fMassJpsi;
336 
337  return *this;
338 }
339 //----------------------------------------------------------------------------
342  if(fV1) { delete fV1; fV1=0; }
343  delete fVertexerTracks;
344  if(fTrackFilter) { delete fTrackFilter; fTrackFilter=0; }
348  if(fCutsD0toKpi) { delete fCutsD0toKpi; fCutsD0toKpi=0; }
349  if(fCutsJpsitoee) { delete fCutsJpsitoee; fCutsJpsitoee=0; }
351  if(fCutsDstoKKpi) { delete fCutsDstoKKpi; fCutsDstoKKpi=0; }
352  if(fCutsLctopKpi) { delete fCutsLctopKpi; fCutsLctopKpi=0; }
353  if(fCutsLctoV0) { delete fCutsLctoV0; fCutsLctoV0=0; }
356  if(fAODMap) { delete [] fAODMap; fAODMap=0; }
357  if(fMassCalc2) { delete fMassCalc2; fMassCalc2=0; }
358  if(fMassCalc3) { delete fMassCalc3; fMassCalc3=0; }
359  if(fMassCalc4) { delete fMassCalc4; fMassCalc4=0; }
360 }
361 //----------------------------------------------------------------------------
364 
365  TList *list = new TList();
366  list->SetOwner();
367  list->SetName("ListOfCuts");
368 
369  if(fCutsD0toKpi) {
371  list->Add(cutsD0toKpi);
372  }
373  if(fCutsJpsitoee) {
375  list->Add(cutsJpsitoee);
376  }
377  if(fCutsDplustoKpipi) {
379  list->Add(cutsDplustoKpipi);
380  }
381  if(fCutsDstoKKpi) {
383  list->Add(cutsDstoKKpi);
384  }
385  if(fCutsLctopKpi) {
387  list->Add(cutsLctopKpi);
388  }
389  if(fCutsLctoV0){
390  AliRDHFCutsLctoV0 *cutsLctoV0 = new AliRDHFCutsLctoV0(*fCutsLctoV0);
391  list->Add(cutsLctoV0);
392  }
393  if(fCutsD0toKpipipi) {
395  list->Add(cutsD0toKpipipi);
396  }
397  if(fCutsDStartoKpipi) {
399  list->Add(cutsDStartoKpipi);
400  }
401 
402  //___ Check consitstency of cuts between vertexer and analysis tasks
403  Bool_t bCutsOk = CheckCutsConsistency();
404  if (bCutsOk == kFALSE) {AliFatal("AliAnalysisVertexingHF::FillListOfCuts vertexing and the analysis task cuts are not consistent!");}
405 
406  // keep a pointer to the list
407  fListOfCuts = list;
408 
409  return list;
410 }
411 //----------------------------------------------------------------------------
413  TClonesArray *aodVerticesHFTClArr,
414  TClonesArray *aodD0toKpiTClArr,
415  TClonesArray *aodJPSItoEleTClArr,
416  TClonesArray *aodCharm3ProngTClArr,
417  TClonesArray *aodCharm4ProngTClArr,
418  TClonesArray *aodDstarTClArr,
419  TClonesArray *aodCascadesTClArr,
420  TClonesArray *aodLikeSign2ProngTClArr,
421  TClonesArray *aodLikeSign3ProngTClArr)
422 {
426  //AliCodeTimerAuto("",0);
427 
428  if(!fMixEvent){
429  TString evtype = event->IsA()->GetName();
430  fInputAOD = ((evtype=="AliAODEvent") ? kTRUE : kFALSE);
431  } // if we do mixing AliVEvent is a AliMixedEvent
432 
433  if(fInputAOD) {
434  AliDebug(2,"Creating HF candidates from AOD");
435  } else {
436  AliDebug(2,"Creating HF candidates from ESD");
437  }
438 
439  if(!aodVerticesHFTClArr) {
440  printf("ERROR: no aodVerticesHFTClArr");
441  return;
442  }
443  if((fD0toKpi || fDstar) && !aodD0toKpiTClArr) {
444  printf("ERROR: no aodD0toKpiTClArr");
445  return;
446  }
447  if(fJPSItoEle && !aodJPSItoEleTClArr) {
448  printf("ERROR: no aodJPSItoEleTClArr");
449  return;
450  }
451  if(f3Prong && !aodCharm3ProngTClArr) {
452  printf("ERROR: no aodCharm3ProngTClArr");
453  return;
454  }
455  if(f4Prong && !aodCharm4ProngTClArr) {
456  printf("ERROR: no aodCharm4ProngTClArr");
457  return;
458  }
459  if(fDstar && !aodDstarTClArr) {
460  printf("ERROR: no aodDstarTClArr");
461  return;
462  }
463  if(fCascades && !aodCascadesTClArr){
464  printf("ERROR: no aodCascadesTClArr ");
465  return;
466  }
467  if(fLikeSign && !aodLikeSign2ProngTClArr) {
468  printf("ERROR: no aodLikeSign2ProngTClArr");
469  return;
470  }
471  if(fLikeSign3prong && f3Prong && !aodLikeSign3ProngTClArr) {
472  printf("ERROR: no aodLikeSign3ProngTClArr");
473  return;
474  }
475 
476  // delete candidates from previous event and create references
477  Int_t iVerticesHF=0,iD0toKpi=0,iJPSItoEle=0,i3Prong=0,i4Prong=0,iDstar=0,iCascades=0,iLikeSign2Prong=0,iLikeSign3Prong=0;
478  aodVerticesHFTClArr->Delete();
479  iVerticesHF = aodVerticesHFTClArr->GetEntriesFast();
480  TClonesArray &verticesHFRef = *aodVerticesHFTClArr;
481  if(fD0toKpi || fDstar) {
482  aodD0toKpiTClArr->Delete();
483  iD0toKpi = aodD0toKpiTClArr->GetEntriesFast();
484  }
485  if(fJPSItoEle) {
486  aodJPSItoEleTClArr->Delete();
487  iJPSItoEle = aodJPSItoEleTClArr->GetEntriesFast();
488  }
489  if(f3Prong) {
490  aodCharm3ProngTClArr->Delete();
491  i3Prong = aodCharm3ProngTClArr->GetEntriesFast();
492  }
493  if(f4Prong) {
494  aodCharm4ProngTClArr->Delete();
495  i4Prong = aodCharm4ProngTClArr->GetEntriesFast();
496  }
497  if(fDstar) {
498  aodDstarTClArr->Delete();
499  iDstar = aodDstarTClArr->GetEntriesFast();
500  }
501  if(fCascades) {
502  aodCascadesTClArr->Delete();
503  iCascades = aodCascadesTClArr->GetEntriesFast();
504  }
505  if(fLikeSign) {
506  aodLikeSign2ProngTClArr->Delete();
507  iLikeSign2Prong = aodLikeSign2ProngTClArr->GetEntriesFast();
508  }
509  if(fLikeSign3prong && f3Prong) {
510  aodLikeSign3ProngTClArr->Delete();
511  iLikeSign3Prong = aodLikeSign3ProngTClArr->GetEntriesFast();
512  }
513 
514  TClonesArray &aodD0toKpiRef = *aodD0toKpiTClArr;
515  TClonesArray &aodJPSItoEleRef = *aodJPSItoEleTClArr;
516  TClonesArray &aodCharm3ProngRef = *aodCharm3ProngTClArr;
517  TClonesArray &aodCharm4ProngRef = *aodCharm4ProngTClArr;
518  TClonesArray &aodDstarRef = *aodDstarTClArr;
519  TClonesArray &aodCascadesRef = *aodCascadesTClArr;
520  TClonesArray &aodLikeSign2ProngRef = *aodLikeSign2ProngTClArr;
521  TClonesArray &aodLikeSign3ProngRef = *aodLikeSign3ProngTClArr;
522 
523 
524  AliAODRecoDecayHF2Prong *io2Prong = 0;
525  AliAODRecoDecayHF3Prong *io3Prong = 0;
526  AliAODRecoDecayHF4Prong *io4Prong = 0;
527  AliAODRecoCascadeHF *ioCascade = 0;
528 
529  Int_t iTrkP1,iTrkP2,iTrkN1,iTrkN2,iTrkSoftPi,trkEntries,iv0,nv0;
530  Double_t xdummy,ydummy,dcap1n1,dcap1n2,dcap2n1,dcap1p2,dcan1n2,dcap2n2,dcaV0,dcaCasc;
531  Bool_t okD0=kFALSE,okJPSI=kFALSE,ok3Prong=kFALSE,ok4Prong=kFALSE;
532  Bool_t okDstar=kFALSE,okD0fromDstar=kFALSE;
533  Bool_t okCascades=kFALSE;
534  AliESDtrack *postrack1 = 0;
535  AliESDtrack *postrack2 = 0;
536  AliESDtrack *negtrack1 = 0;
537  AliESDtrack *negtrack2 = 0;
538  AliESDtrack *trackPi = 0;
539  Double_t mompos1[3],mompos2[3],momneg1[3],momneg2[3];
540  // AliESDtrack *posV0track = 0;
541  // AliESDtrack *negV0track = 0;
542  Float_t dcaMax = fCutsD0toKpi->GetDCACut();
543  if(fCutsJpsitoee) dcaMax=TMath::Max(dcaMax,fCutsJpsitoee->GetDCACut());
544  if(fCutsDplustoKpipi) dcaMax=TMath::Max(dcaMax,fCutsDplustoKpipi->GetDCACut());
545  if(fCutsDstoKKpi) dcaMax=TMath::Max(dcaMax,fCutsDstoKKpi->GetDCACut());
546  if(fCutsLctopKpi) dcaMax=TMath::Max(dcaMax,fCutsLctopKpi->GetDCACut());
547  if(fCutsD0toKpipipi) dcaMax=TMath::Max(dcaMax,fCutsD0toKpipipi->GetDCACut());
548  if(fCutsDStartoKpipi) dcaMax=TMath::Max(dcaMax,fCutsDStartoKpipi->GetDCACut());
549 
550  AliDebug(2,Form(" dca cut set to %f cm",dcaMax));
551 
552 
553  // get Bz
554  fBzkG = (Double_t)event->GetMagneticField();
555  if(!fVertexerTracks){
556  fVertexerTracks=new AliVertexerTracks(fBzkG);
557  }else{
558  Double_t oldField=fVertexerTracks->GetFieldkG();
559  if(oldField!=fBzkG) fVertexerTracks->SetFieldkG(fBzkG);
560  }
561 
562  trkEntries = (Int_t)event->GetNumberOfTracks();
563  AliDebug(1,Form(" Number of tracks: %d",trkEntries));
564  fnTrksTotal += trkEntries;
565 
566  nv0 = (Int_t)event->GetNumberOfV0s();
567  AliDebug(1,Form(" Number of V0s: %d",nv0));
568 
569  if( trkEntries<2 && (trkEntries<1 || nv0<1) ) {
570  AliDebug(1,Form(" Not enough tracks: %d",trkEntries));
571  return;
572  }
573 
574  // event selection + PID configuration
575  if(!fCutsD0toKpi->IsEventSelected(event)) return;
580  if(fCutsLctoV0) fCutsLctoV0->SetupPID(event);
583 
584  // call function that applies sigle-track selection,
585  // for displaced tracks and soft pions (both charges) for D*,
586  // and retrieves primary vertex
587  TObjArray seleTrksArray(trkEntries);
588  TObjArray tracksAtVertex(trkEntries);
589  UChar_t *seleFlags = new UChar_t[trkEntries]; // bit 0: displaced, bit 1: softpi, bit 2: 3 prong, bits 3-4-5: for PID
590  Int_t nSeleTrks=0;
591  Int_t *evtNumber = new Int_t[trkEntries];
592  SelectTracksAndCopyVertex(event,trkEntries,seleTrksArray,tracksAtVertex,nSeleTrks,seleFlags,evtNumber);
593 
594  AliDebug(1,Form(" Selected tracks: %d",nSeleTrks));
595  fnSeleTrksTotal += nSeleTrks;
596 
597 
598  TObjArray *twoTrackArray1 = new TObjArray(2);
599  TObjArray *twoTrackArray2 = new TObjArray(2);
600  TObjArray *twoTrackArrayV0 = new TObjArray(2);
601  TObjArray *twoTrackArrayCasc = new TObjArray(2);
602  TObjArray *threeTrackArray = new TObjArray(3);
603  TObjArray *fourTrackArray = new TObjArray(4);
604 
605  Double_t dispersion;
606  Bool_t isLikeSign2Prong=kFALSE,isLikeSign3Prong=kFALSE;
607 
608  AliAODRecoDecayHF *rd = 0;
609  AliAODRecoCascadeHF *rc = 0;
610  AliAODv0 *v0 = 0;
611  AliESDv0 *esdV0 = 0;
612 
613  Bool_t massCutOK=kTRUE;
614 
615  // LOOP ON POSITIVE TRACKS
616  for(iTrkP1=0; iTrkP1<nSeleTrks; iTrkP1++) {
617 
618  //if(iTrkP1%1==0) AliDebug(1,Form(" 1st loop on pos: track number %d of %d",iTrkP1,nSeleTrks));
619  //if(iTrkP1%1==0) printf(" 1st loop on pos: track number %d of %d\n",iTrkP1,nSeleTrks);
620 
621  // get track from tracks array
622  postrack1 = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkP1);
623  if(!TESTBIT(seleFlags[iTrkP1],kBitDispl)) continue;
624  postrack1->GetPxPyPz(mompos1);
625 
626  // Make cascades with V0+track
627  //
628  if(fCascades) {
629  // loop on V0's
630  for(iv0=0; iv0<nv0; iv0++){
631 
632  //AliDebug(1,Form(" loop on v0s for track number %d and v0 number %d",iTrkP1,iv0));
633 
634  if ( fUsePIDforLc2V0 && !TESTBIT(seleFlags[iTrkP1],kBitProtonCompat) ) continue; //clm
635 
636  // Get the V0
637  if(fInputAOD) {
638  v0 = ((AliAODEvent*)event)->GetV0(iv0);
639  } else {
640  esdV0 = ((AliESDEvent*)event)->GetV0(iv0);
641  }
642  if ( (!v0 || !v0->IsA()->InheritsFrom("AliAODv0") ) &&
643  (!esdV0 || !esdV0->IsA()->InheritsFrom("AliESDv0") ) ) continue;
644 
645  if ( v0 && ((v0->GetOnFlyStatus() == kTRUE && fV0TypeForCascadeVertex == AliRDHFCuts::kOnlyOfflineV0s) ||
646  (v0->GetOnFlyStatus() == kFALSE && fV0TypeForCascadeVertex == AliRDHFCuts::kOnlyOnTheFlyV0s)) ) continue;
647 
648  if ( esdV0 && ((esdV0->GetOnFlyStatus() == kTRUE && fV0TypeForCascadeVertex == AliRDHFCuts::kOnlyOfflineV0s) ||
649  ( esdV0->GetOnFlyStatus() == kFALSE && fV0TypeForCascadeVertex == AliRDHFCuts::kOnlyOnTheFlyV0s)) ) continue;
650 
651  // Get the tracks that form the V0
652  // ( parameters at primary vertex )
653  // and define an AliExternalTrackParam out of them
654  AliExternalTrackParam * posV0track;
655  AliExternalTrackParam * negV0track;
656 
657  if(fInputAOD){
658  AliAODTrack *posVV0track = (AliAODTrack*)(v0->GetDaughter(0));
659  AliAODTrack *negVV0track = (AliAODTrack*)(v0->GetDaughter(1));
660  if( !posVV0track || !negVV0track ) continue;
661  //
662  // Apply some basic V0 daughter criteria
663  //
664  // bachelor must not be a v0-track
665  if (posVV0track->GetID() == postrack1->GetID() ||
666  negVV0track->GetID() == postrack1->GetID()) continue;
667  // reject like-sign v0
668  if ( posVV0track->Charge() == negVV0track->Charge() ) continue;
669  // avoid ghost TPC tracks
670  if(!(posVV0track->GetStatus() & AliESDtrack::kTPCrefit) ||
671  !(negVV0track->GetStatus() & AliESDtrack::kTPCrefit)) continue;
672  // Get AliExternalTrackParam out of the AliAODTracks
673  Double_t xyz[3], pxpypz[3], cv[21]; Short_t sign;
674  posVV0track->PxPyPz(pxpypz); posVV0track->XvYvZv(xyz);
675  posVV0track->GetCovarianceXYZPxPyPz(cv); sign=posVV0track->Charge();
676  posV0track = new AliExternalTrackParam(xyz,pxpypz,cv,sign);
677  negVV0track->PxPyPz(pxpypz); negVV0track->XvYvZv(xyz);
678  negVV0track->GetCovarianceXYZPxPyPz(cv); sign=negVV0track->Charge();
679  negV0track = new AliExternalTrackParam(xyz,pxpypz,cv,sign);
680  } else {
681  AliESDtrack *posVV0track = (AliESDtrack*)(event->GetTrack( esdV0->GetPindex() ));
682  AliESDtrack *negVV0track = (AliESDtrack*)(event->GetTrack( esdV0->GetNindex() ));
683  if( !posVV0track || !negVV0track ) continue;
684  //
685  // Apply some basic V0 daughter criteria
686  //
687  // bachelor must not be a v0-track
688  if (posVV0track->GetID() == postrack1->GetID() ||
689  negVV0track->GetID() == postrack1->GetID()) continue;
690  // reject like-sign v0
691  if ( posVV0track->Charge() == negVV0track->Charge() ) continue;
692  // avoid ghost TPC tracks
693  if(!(posVV0track->GetStatus() & AliESDtrack::kTPCrefit) ||
694  !(negVV0track->GetStatus() & AliESDtrack::kTPCrefit)) continue;
695  // reject kinks (only necessary on AliESDtracks)
696  if (posVV0track->GetKinkIndex(0)>0 || negVV0track->GetKinkIndex(0)>0) continue;
697  // Get AliExternalTrackParam out of the AliESDtracks
698  posV0track = new AliExternalTrackParam(*posVV0track);
699  negV0track = new AliExternalTrackParam(*negVV0track);
700 
701  // Define the AODv0 from ESDv0 if reading ESDs
702  v0 = TransformESDv0toAODv0(esdV0,twoTrackArrayV0);
703  }
704  if( !posV0track || !negV0track ){
705  AliDebug(1,Form(" Couldn't get the V0 daughters"));
706  continue;
707  }
708 
709  // fill in the v0 two-external-track-param array
710  twoTrackArrayV0->AddAt(posV0track,0);
711  twoTrackArrayV0->AddAt(negV0track,1);
712 
713  // Get the V0 dca
714  dcaV0 = v0->DcaV0Daughters();
715 
716  // Define the V0 (neutral) track
717  AliNeutralTrackParam *trackV0=NULL;
718  if(fInputAOD) {
719  const AliVTrack *trackVV0 = dynamic_cast<const AliVTrack*>(v0);
720  if(trackVV0) trackV0 = new AliNeutralTrackParam(trackVV0);
721  } else {
722  Double_t xyz[3], pxpypz[3];
723  esdV0->XvYvZv(xyz);
724  esdV0->PxPyPz(pxpypz);
725  Double_t cv[21]; for(int i=0; i<21; i++) cv[i]=0;
726  trackV0 = new AliNeutralTrackParam(xyz,pxpypz,cv,0);
727  }
728 
729 
730  // Fill in the object array to create the cascade
731  twoTrackArrayCasc->AddAt(postrack1,0);
732  twoTrackArrayCasc->AddAt(trackV0,1);
733  // Compute the cascade vertex
734  AliAODVertex *vertexCasc = 0;
736  // DCA between the two tracks
737  dcaCasc = postrack1->GetDCA(trackV0,fBzkG,xdummy,ydummy);
738  // Vertexing+
739  vertexCasc = ReconstructSecondaryVertex(twoTrackArrayCasc,dispersion,kFALSE);
740  } else {
741  // assume Cascade decays at the primary vertex
742  Double_t pos[3],cov[6],chi2perNDF;
743  fV1->GetXYZ(pos);
744  fV1->GetCovMatrix(cov);
745  chi2perNDF = fV1->GetChi2toNDF();
746  vertexCasc = new AliAODVertex(pos,cov,chi2perNDF,0x0,-1,AliAODVertex::kUndef,2);
747  dcaCasc = 0.;
748  }
749  if(!vertexCasc) {
750  delete posV0track; posV0track=NULL;
751  delete negV0track; negV0track=NULL;
752  delete trackV0; trackV0=NULL;
753  if(!fInputAOD) {delete v0; v0=NULL;}
754  twoTrackArrayV0->Clear();
755  twoTrackArrayCasc->Clear();
756  continue;
757  }
758 
759  // Create and store the Cascade if passed the cuts
760  ioCascade = MakeCascade(twoTrackArrayCasc,event,vertexCasc,v0,dcaCasc,okCascades);
761  if(okCascades && ioCascade) {
762  //AliDebug(1,Form("Storing a cascade object... "));
763  // add the vertex and the cascade to the AOD
764  rc = new(aodCascadesRef[iCascades++])AliAODRecoCascadeHF(*ioCascade);
765  if(fMakeReducedRHF){
766  UShort_t id[2]={(UShort_t)postrack1->GetID(),(UShort_t)iv0};
767  rc->SetProngIDs(2,id);
768  rc->DeleteRecoD();
769  }else{
770  AliAODVertex *vCasc = new(verticesHFRef[iVerticesHF++])AliAODVertex(*vertexCasc);
771  rc->SetSecondaryVtx(vCasc);
772  vCasc->SetParent(rc);
773  if(!fInputAOD) vCasc->AddDaughter(v0); // just to fill ref #0 ??
774  AddRefs(vCasc,rc,event,twoTrackArrayCasc); // add the track (proton)
775  vCasc->AddDaughter(v0); // fill the 2prong V0
776  }
777  rc->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
778  }
779 
780 
781  // Clean up
782  delete posV0track; posV0track=NULL;
783  delete negV0track; negV0track=NULL;
784  delete trackV0; trackV0=NULL;
785  twoTrackArrayV0->Clear();
786  twoTrackArrayCasc->Clear();
787  if(ioCascade) { delete ioCascade; ioCascade=NULL; }
788  if(vertexCasc) { delete vertexCasc; vertexCasc=NULL; }
789  if(!fInputAOD) {delete v0; v0=NULL;}
790 
791  } // end loop on V0's
792  }
793 
794  // If there is less than 2 particles continue
795  if(trkEntries<2) {
796  AliDebug(1,Form(" Not enough tracks: %d",trkEntries));
797  continue;
798  }
799 
800  if(postrack1->Charge()<0 && !fLikeSign) continue;
801 
802  // LOOP ON NEGATIVE TRACKS
803  for(iTrkN1=0; iTrkN1<nSeleTrks; iTrkN1++) {
804 
805  //if(iTrkN1%1==0) AliDebug(1,Form(" 1st loop on neg: track number %d of %d",iTrkN1,nSeleTrks));
806  //if(iTrkN1%1==0) printf(" 1st loop on neg: track number %d of %d\n",iTrkN1,nSeleTrks);
807 
808  if(iTrkN1==iTrkP1) continue;
809 
810  // get track from tracks array
811  negtrack1 = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkN1);
812 
813  if(negtrack1->Charge()>0 && !fLikeSign) continue;
814 
815  if(!TESTBIT(seleFlags[iTrkN1],kBitDispl)) continue;
816 
817  if(fMixEvent) {
818  if(evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
819  }
820 
821  if(postrack1->Charge()==negtrack1->Charge()) { // like-sign
822  isLikeSign2Prong=kTRUE;
823  if(!fLikeSign) continue;
824  if(iTrkN1<iTrkP1) continue; // this is needed to avoid double-counting of like-sign
825  } else { // unlike-sign
826  isLikeSign2Prong=kFALSE;
827  if(postrack1->Charge()<0 || negtrack1->Charge()>0) continue; // this is needed to avoid double-counting of unlike-sign
828  if(fMixEvent) {
829  if(evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
830  }
831 
832  }
833 
834  // back to primary vertex
835  // postrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
836  // negtrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
837  SetParametersAtVertex(postrack1,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkP1));
838  SetParametersAtVertex(negtrack1,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkN1));
839  negtrack1->GetPxPyPz(momneg1);
840 
841  // DCA between the two tracks
842  dcap1n1 = postrack1->GetDCA(negtrack1,fBzkG,xdummy,ydummy);
843  if(dcap1n1>dcaMax) { negtrack1=0; continue; }
844 
845  // Vertexing
846  twoTrackArray1->AddAt(postrack1,0);
847  twoTrackArray1->AddAt(negtrack1,1);
848  AliAODVertex *vertexp1n1 = ReconstructSecondaryVertex(twoTrackArray1,dispersion);
849  if(!vertexp1n1) {
850  twoTrackArray1->Clear();
851  negtrack1=0;
852  continue;
853  }
854  // 2 prong candidate
855  if(fD0toKpi || fJPSItoEle || fDstar || fLikeSign) {
856 
857  io2Prong = Make2Prong(twoTrackArray1,event,vertexp1n1,dcap1n1,okD0,okJPSI,okD0fromDstar);
858 
859  if((fD0toKpi && okD0) || (fJPSItoEle && okJPSI) || (isLikeSign2Prong && (okD0 || okJPSI))) {
860  // add the vertex and the decay to the AOD
861  AliAODVertex *v2Prong =0x0;
862  if(!fMakeReducedRHF)v2Prong = new(verticesHFRef[iVerticesHF++])AliAODVertex(*vertexp1n1);
863  if(!isLikeSign2Prong) {
864  if(okD0) {
865  rd = new(aodD0toKpiRef[iD0toKpi++])AliAODRecoDecayHF2Prong(*io2Prong);
867 
868  if(fMakeReducedRHF){
869  rd->DeleteRecoD();
870  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
871  }else{
872  rd->SetSecondaryVtx(v2Prong);
873  v2Prong->SetParent(rd);
874  AddRefs(v2Prong,rd,event,twoTrackArray1);
875  }
876  }
877  if(okJPSI) {
878  rd = new(aodJPSItoEleRef[iJPSItoEle++])AliAODRecoDecayHF2Prong(*io2Prong);
879  if(fMakeReducedRHF){
880  rd->DeleteRecoD();
881  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
882  }else{
883  if(!okD0) v2Prong->SetParent(rd); // it cannot have two mothers ...
884  AddRefs(v2Prong,rd,event,twoTrackArray1);
885  }
886  }
887  } else { // isLikeSign2Prong
888  rd = new(aodLikeSign2ProngRef[iLikeSign2Prong++])AliAODRecoDecayHF2Prong(*io2Prong);
889  //Set selection bit for PID
891  if(fMakeReducedRHF){
892  rd->DeleteRecoD();
893  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
894  }else{
895  rd->SetSecondaryVtx(v2Prong);
896  v2Prong->SetParent(rd);
897  AddRefs(v2Prong,rd,event,twoTrackArray1);
898  }
899  }
900  }
901  // D* candidates
902  if(fDstar && okD0fromDstar && !isLikeSign2Prong) {
903  // write references in io2Prong
904  if(fInputAOD) {
905  AddDaughterRefs(vertexp1n1,event,twoTrackArray1);
906  } else {
907  vertexp1n1->AddDaughter(postrack1);
908  vertexp1n1->AddDaughter(negtrack1);
909  }
910  io2Prong->SetSecondaryVtx(vertexp1n1);
911  //printf("---> %d %d %d %d %d\n",vertexp1n1->GetNDaughters(),iTrkP1,iTrkN1,postrack1->Charge(),negtrack1->Charge());
912  // create a track from the D0
913  AliNeutralTrackParam *trackD0 = new AliNeutralTrackParam(io2Prong);
914 
915  // LOOP ON TRACKS THAT PASSED THE SOFT PION CUTS
916  for(iTrkSoftPi=0; iTrkSoftPi<nSeleTrks; iTrkSoftPi++) {
917 
918  if(iTrkSoftPi==iTrkP1 || iTrkSoftPi==iTrkN1) continue;
919 
920  if(!TESTBIT(seleFlags[iTrkSoftPi],kBitSoftPi)) continue;
921 
922  if(fMixEvent) {
923  if(evtNumber[iTrkP1]==evtNumber[iTrkSoftPi] ||
924  evtNumber[iTrkN1]==evtNumber[iTrkSoftPi] ||
925  evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
926  }
927 
928  //if(iTrkSoftPi%1==0) AliDebug(1,Form(" 1st loop on pi_s: track number %d of %d",iTrkSoftPi,nSeleTrks));
929 
930  trackD0->PropagateToDCA(fV1,fBzkG,kVeryBig);
931  if(trackD0->GetSigmaY2()<0. || trackD0->GetSigmaZ2()<0.) continue; // this is insipired by the AliITStrackV2::Invariant() checks
932 
933  // get track from tracks array
934  trackPi = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkSoftPi);
935  // trackPi->PropagateToDCA(fV1,fBzkG,kVeryBig);
936  SetParametersAtVertex(trackPi,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkSoftPi));
937  twoTrackArrayCasc->AddAt(trackPi,0);
938  twoTrackArrayCasc->AddAt(trackD0,1);
939  if(!SelectInvMassAndPtDstarD0pi(twoTrackArrayCasc)){
940  twoTrackArrayCasc->Clear();
941  trackPi=0;
942  continue;
943  }
944 
945  AliAODVertex *vertexCasc = 0;
946 
947  if(fFindVertexForDstar) {
948  // DCA between the two tracks
949  dcaCasc = trackPi->GetDCA(trackD0,fBzkG,xdummy,ydummy);
950  // Vertexing
951  vertexCasc = ReconstructSecondaryVertex(twoTrackArrayCasc,dispersion,kFALSE);
952  } else {
953  // assume Dstar decays at the primary vertex
954  Double_t pos[3],cov[6],chi2perNDF;
955  fV1->GetXYZ(pos);
956  fV1->GetCovMatrix(cov);
957  chi2perNDF = fV1->GetChi2toNDF();
958  vertexCasc = new AliAODVertex(pos,cov,chi2perNDF,0x0,-1,AliAODVertex::kUndef,2);
959  dcaCasc = 0.;
960  }
961  if(!vertexCasc) {
962  twoTrackArrayCasc->Clear();
963  trackPi=0;
964  continue;
965  }
966 
967  ioCascade = MakeCascade(twoTrackArrayCasc,event,vertexCasc,io2Prong,dcaCasc,okDstar);
968  if(okDstar) {
969  // add the D0 to the AOD (if not already done)
970  if(!okD0) {
971  rd = new(aodD0toKpiRef[iD0toKpi++])AliAODRecoDecayHF2Prong(*io2Prong);
972  if(fMakeReducedRHF){
973  rd->DeleteRecoD();
974  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
975  }else{
976  AliAODVertex *v2Prong = new (verticesHFRef[iVerticesHF++])AliAODVertex(*vertexp1n1);
977  rd->SetSecondaryVtx(v2Prong);
978  v2Prong->SetParent(rd);
979  AddRefs(v2Prong,rd,event,twoTrackArray1);
980  }
981  okD0=kTRUE; // this is done to add it only once
982  }
983  // add the vertex and the cascade to the AOD
984  rc = new(aodDstarRef[iDstar++])AliAODRecoCascadeHF(*ioCascade);
985  // Set selection bit for PID
987  AliAODVertex *vCasc = 0x0;
988  if(fMakeReducedRHF){
989  //assign a ID to the D0 candidate, daughter of the Cascade. ID = position in the D0toKpi array
990  UShort_t idCasc[2]={(UShort_t)trackPi->GetID(),(UShort_t)(iD0toKpi-1)};
991  rc->SetProngIDs(2,idCasc);
992  rc->DeleteRecoD();
993  rc->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
994  }else{
995  AliAODVertex *vCasc = new(verticesHFRef[iVerticesHF++])AliAODVertex(*vertexCasc);
996  rc->SetSecondaryVtx(vCasc);
997  vCasc->SetParent(rc);
998  if(!fInputAOD) vCasc->AddDaughter(rd); // just to fill ref #0
999  AddRefs(vCasc,rc,event,twoTrackArrayCasc);
1000  vCasc->AddDaughter(rd); // add the D0 (in ref #1)
1001  }
1002  }
1003  twoTrackArrayCasc->Clear();
1004  trackPi=0;
1005  if(ioCascade) {delete ioCascade; ioCascade=NULL;}
1006  delete vertexCasc; vertexCasc=NULL;
1007  } // end loop on soft pi tracks
1008 
1009  if(trackD0) {delete trackD0; trackD0=NULL;}
1010 
1011  }
1012  if(io2Prong) {delete io2Prong; io2Prong=NULL;}
1013  }
1014 
1015  twoTrackArray1->Clear();
1016  if( (!f3Prong && !f4Prong) ||
1017  (isLikeSign2Prong && !f3Prong) ) {
1018  negtrack1=0;
1019  delete vertexp1n1;
1020  continue;
1021  }
1022 
1023 
1024  // 2nd LOOP ON POSITIVE TRACKS
1025  for(iTrkP2=iTrkP1+1; iTrkP2<nSeleTrks; iTrkP2++) {
1026 
1027  if(iTrkP2==iTrkP1 || iTrkP2==iTrkN1) continue;
1028 
1029  //if(iTrkP2%1==0) AliDebug(1,Form(" 2nd loop on pos: track number %d of %d",iTrkP2,nSeleTrks));
1030 
1031  // get track from tracks array
1032  postrack2 = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkP2);
1033 
1034  if(postrack2->Charge()<0) continue;
1035 
1036  if(!TESTBIT(seleFlags[iTrkP2],kBitDispl)) continue;
1037 
1038  // Check single tracks cuts specific for 3 prongs
1039  if(!TESTBIT(seleFlags[iTrkP2],kBit3Prong)) continue;
1040  if(!TESTBIT(seleFlags[iTrkP1],kBit3Prong)) continue;
1041  if(!TESTBIT(seleFlags[iTrkN1],kBit3Prong)) continue;
1042 
1043  if(fMixEvent) {
1044  if(evtNumber[iTrkP1]==evtNumber[iTrkP2] ||
1045  evtNumber[iTrkN1]==evtNumber[iTrkP2] ||
1046  evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
1047  }
1048 
1049  if(isLikeSign2Prong) { // like-sign pair -> have to build only like-sign triplet
1050  if(!fLikeSign3prong) continue;
1051  if(postrack1->Charge()>0) { // ok: like-sign triplet (+++)
1052  isLikeSign3Prong=kTRUE;
1053  } else { // not ok
1054  continue;
1055  }
1056  } else { // normal triplet (+-+)
1057  isLikeSign3Prong=kFALSE;
1058  if(fMixEvent) {
1059  if(evtNumber[iTrkP1]==evtNumber[iTrkP2] ||
1060  evtNumber[iTrkN1]==evtNumber[iTrkP2] ||
1061  evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
1062  }
1063  }
1064 
1066  if(!TESTBIT(seleFlags[iTrkN1],kBitKaonCompat)) continue;
1067  }
1068  Bool_t okForLcTopKpi=kTRUE;
1069  Int_t pidLcStatus=3; // 3= OK as pKpi and Kpipi
1070  if(fUsePIDforLc>0){
1071  if(!TESTBIT(seleFlags[iTrkP1],kBitProtonCompat) &&
1072  !TESTBIT(seleFlags[iTrkP2],kBitProtonCompat) ){
1073  okForLcTopKpi=kFALSE;
1074  pidLcStatus=0;
1075  }
1076  if(okForLcTopKpi && fUsePIDforLc>1){
1077  okForLcTopKpi=kFALSE;
1078  pidLcStatus=0;
1079  if(TESTBIT(seleFlags[iTrkP1],kBitProtonCompat) &&
1080  TESTBIT(seleFlags[iTrkP2],kBitPionCompat) ){
1081  okForLcTopKpi=kTRUE;
1082  pidLcStatus+=1; // 1= OK as pKpi
1083  }
1084  if(TESTBIT(seleFlags[iTrkP2],kBitProtonCompat) &&
1085  TESTBIT(seleFlags[iTrkP1],kBitPionCompat) ){
1086  okForLcTopKpi=kTRUE;
1087  pidLcStatus+=2; // 2= OK as piKp
1088  }
1089  }
1090  }
1091  Bool_t okForDsToKKpi=kTRUE;
1092  if(fUseKaonPIDforDs){
1093  if(!TESTBIT(seleFlags[iTrkP1],kBitKaonCompat) &&
1094  !TESTBIT(seleFlags[iTrkP2],kBitKaonCompat) ) okForDsToKKpi=kFALSE;
1095  }
1096  // back to primary vertex
1097  // postrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
1098  // postrack2->PropagateToDCA(fV1,fBzkG,kVeryBig);
1099  // negtrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
1100  SetParametersAtVertex(postrack1,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkP1));
1101  SetParametersAtVertex(negtrack1,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkN1));
1102  SetParametersAtVertex(postrack2,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkP2));
1103 
1104  //printf("********** %d %d %d\n",postrack1->GetID(),postrack2->GetID(),negtrack1->GetID());
1105 
1106  dcap2n1 = postrack2->GetDCA(negtrack1,fBzkG,xdummy,ydummy);
1107  if(dcap2n1>dcaMax) { postrack2=0; continue; }
1108  dcap1p2 = postrack2->GetDCA(postrack1,fBzkG,xdummy,ydummy);
1109  if(dcap1p2>dcaMax) { postrack2=0; continue; }
1110 
1111  // check invariant mass cuts for D+,Ds,Lc
1112  massCutOK=kTRUE;
1113  if(f3Prong) {
1114  if(postrack2->Charge()>0) {
1115  threeTrackArray->AddAt(postrack1,0);
1116  threeTrackArray->AddAt(negtrack1,1);
1117  threeTrackArray->AddAt(postrack2,2);
1118  } else {
1119  threeTrackArray->AddAt(negtrack1,0);
1120  threeTrackArray->AddAt(postrack1,1);
1121  threeTrackArray->AddAt(postrack2,2);
1122  }
1124  postrack2->GetPxPyPz(mompos2);
1125  Double_t pxDau[3]={mompos1[0],momneg1[0],mompos2[0]};
1126  Double_t pyDau[3]={mompos1[1],momneg1[1],mompos2[1]};
1127  Double_t pzDau[3]={mompos1[2],momneg1[2],mompos2[2]};
1128  // massCutOK = SelectInvMassAndPt3prong(threeTrackArray);
1129  massCutOK = SelectInvMassAndPt3prong(pxDau,pyDau,pzDau,pidLcStatus);
1130  }
1131  }
1132 
1133  if(f3Prong && !massCutOK) {
1134  threeTrackArray->Clear();
1135  if(!f4Prong) {
1136  postrack2=0;
1137  continue;
1138  }
1139  }
1140 
1141  // Vertexing
1142  twoTrackArray2->AddAt(postrack2,0);
1143  twoTrackArray2->AddAt(negtrack1,1);
1144  AliAODVertex *vertexp2n1 = ReconstructSecondaryVertex(twoTrackArray2,dispersion);
1145  if(!vertexp2n1) {
1146  twoTrackArray2->Clear();
1147  postrack2=0;
1148  continue;
1149  }
1150 
1151  // 3 prong candidates
1152  if(f3Prong && massCutOK) {
1153 
1154  AliAODVertex* secVert3PrAOD = ReconstructSecondaryVertex(threeTrackArray,dispersion);
1155  io3Prong = Make3Prong(threeTrackArray,event,secVert3PrAOD,dispersion,vertexp1n1,vertexp2n1,dcap1n1,dcap2n1,dcap1p2,okForLcTopKpi,okForDsToKKpi,ok3Prong);
1156  if(ok3Prong) {
1157  AliAODVertex *v3Prong=0x0;
1158  if(!fMakeReducedRHF)v3Prong = new (verticesHFRef[iVerticesHF++])AliAODVertex(*secVert3PrAOD);
1159  if(!isLikeSign3Prong) {
1160  rd = new(aodCharm3ProngRef[i3Prong++])AliAODRecoDecayHF3Prong(*io3Prong);
1161  // Set selection bit for PID
1165  if(fMakeReducedRHF){
1166  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
1167  ((AliAODRecoDecayHF3Prong*)rd)->DeleteRecoD();
1168  }else{
1169  v3Prong = new (verticesHFRef[iVerticesHF++])AliAODVertex(*secVert3PrAOD);
1170  rd->SetSecondaryVtx(v3Prong);
1171  v3Prong->SetParent(rd);
1172  AddRefs(v3Prong,rd,event,threeTrackArray);
1173  }
1174  } else { // isLikeSign3Prong
1175  if(fLikeSign3prong){
1176  rd = new(aodLikeSign3ProngRef[iLikeSign3Prong++])AliAODRecoDecayHF3Prong(*io3Prong);
1177  // Set selection bit for PID
1181  if(fMakeReducedRHF){
1182  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
1183  ((AliAODRecoDecayHF3Prong*)rd)->DeleteRecoD();
1184  }else{
1185  rd->SetSecondaryVtx(v3Prong);
1186  v3Prong->SetParent(rd);
1187  AddRefs(v3Prong,rd,event,threeTrackArray);
1188  }
1189  }
1190  }
1191 
1192  }
1193  if(io3Prong) {delete io3Prong; io3Prong=NULL;}
1194  if(secVert3PrAOD) {delete secVert3PrAOD; secVert3PrAOD=NULL;}
1195  }
1196 
1197  // 4 prong candidates
1198  if(f4Prong
1199  // don't make 4 prong with like-sign pairs and triplets
1200  && !isLikeSign2Prong && !isLikeSign3Prong
1201  // track-to-track dca cuts already now
1202  && dcap1n1 < fCutsD0toKpipipi->GetDCACut()
1203  && dcap2n1 < fCutsD0toKpipipi->GetDCACut()) {
1204  // back to primary vertex
1205  // postrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
1206  // postrack2->PropagateToDCA(fV1,fBzkG,kVeryBig);
1207  // negtrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
1208  SetParametersAtVertex(postrack1,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkP1));
1209  SetParametersAtVertex(negtrack1,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkN1));
1210  SetParametersAtVertex(postrack2,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkP2));
1211 
1212  // Vertexing for these 3 (can be taken from above?)
1213  threeTrackArray->AddAt(postrack1,0);
1214  threeTrackArray->AddAt(negtrack1,1);
1215  threeTrackArray->AddAt(postrack2,2);
1216  AliAODVertex* vertexp1n1p2 = ReconstructSecondaryVertex(threeTrackArray,dispersion);
1217 
1218  // 3rd LOOP ON NEGATIVE TRACKS (for 4 prong)
1219  for(iTrkN2=iTrkN1+1; iTrkN2<nSeleTrks; iTrkN2++) {
1220 
1221  if(iTrkN2==iTrkP1 || iTrkN2==iTrkP2 || iTrkN2==iTrkN1) continue;
1222 
1223  //if(iTrkN2%1==0) AliDebug(1,Form(" 3rd loop on neg: track number %d of %d",iTrkN2,nSeleTrks));
1224 
1225  // get track from tracks array
1226  negtrack2 = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkN2);
1227 
1228  if(negtrack2->Charge()>0) continue;
1229 
1230  if(!TESTBIT(seleFlags[iTrkN2],kBitDispl)) continue;
1231  if(fMixEvent){
1232  if(evtNumber[iTrkP1]==evtNumber[iTrkN2] ||
1233  evtNumber[iTrkN1]==evtNumber[iTrkN2] ||
1234  evtNumber[iTrkP2]==evtNumber[iTrkN2] ||
1235  evtNumber[iTrkP1]==evtNumber[iTrkN1] ||
1236  evtNumber[iTrkP1]==evtNumber[iTrkP2] ||
1237  evtNumber[iTrkN1]==evtNumber[iTrkP2]) continue;
1238  }
1239 
1240  // back to primary vertex
1241  // postrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
1242  // postrack2->PropagateToDCA(fV1,fBzkG,kVeryBig);
1243  // negtrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
1244  // negtrack2->PropagateToDCA(fV1,fBzkG,kVeryBig);
1245  SetParametersAtVertex(postrack1,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkP1));
1246  SetParametersAtVertex(negtrack1,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkN1));
1247  SetParametersAtVertex(postrack2,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkP2));
1248  SetParametersAtVertex(negtrack2,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkN2));
1249 
1250  dcap1n2 = postrack1->GetDCA(negtrack2,fBzkG,xdummy,ydummy);
1251  if(dcap1n2 > fCutsD0toKpipipi->GetDCACut()) { negtrack2=0; continue; }
1252  dcap2n2 = postrack2->GetDCA(negtrack2,fBzkG,xdummy,ydummy);
1253  if(dcap2n2 > fCutsD0toKpipipi->GetDCACut()) { negtrack2=0; continue; }
1254 
1255 
1256  fourTrackArray->AddAt(postrack1,0);
1257  fourTrackArray->AddAt(negtrack1,1);
1258  fourTrackArray->AddAt(postrack2,2);
1259  fourTrackArray->AddAt(negtrack2,3);
1260 
1261  // check invariant mass cuts for D0
1262  massCutOK=kTRUE;
1264  massCutOK = SelectInvMassAndPt4prong(fourTrackArray);
1265 
1266  if(!massCutOK) {
1267  fourTrackArray->Clear();
1268  negtrack2=0;
1269  continue;
1270  }
1271 
1272  // Vertexing
1273  AliAODVertex* secVert4PrAOD = ReconstructSecondaryVertex(fourTrackArray,dispersion);
1274  io4Prong = Make4Prong(fourTrackArray,event,secVert4PrAOD,vertexp1n1,vertexp1n1p2,dcap1n1,dcap1n2,dcap2n1,dcap2n2,ok4Prong);
1275  if(ok4Prong) {
1276  rd = new(aodCharm4ProngRef[i4Prong++])AliAODRecoDecayHF4Prong(*io4Prong);
1277  if(fMakeReducedRHF){
1278  rd->DeleteRecoD();
1279  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
1280  }else{
1281  AliAODVertex *v4Prong = new(verticesHFRef[iVerticesHF++])AliAODVertex(*secVert4PrAOD);
1282  rd->SetSecondaryVtx(v4Prong);
1283  v4Prong->SetParent(rd);
1284  AddRefs(v4Prong,rd,event,fourTrackArray);
1285  }
1286  }
1287 
1288  if(io4Prong) {delete io4Prong; io4Prong=NULL;}
1289  if(secVert4PrAOD) {delete secVert4PrAOD; secVert4PrAOD=NULL;}
1290  fourTrackArray->Clear();
1291  negtrack2 = 0;
1292 
1293  } // end loop on negative tracks
1294 
1295  threeTrackArray->Clear();
1296  delete vertexp1n1p2;
1297 
1298  }
1299 
1300  postrack2 = 0;
1301  delete vertexp2n1;
1302 
1303  } // end 2nd loop on positive tracks
1304 
1305  twoTrackArray2->Clear();
1306 
1307  // 2nd LOOP ON NEGATIVE TRACKS (for 3 prong -+-)
1308  for(iTrkN2=iTrkN1+1; iTrkN2<nSeleTrks; iTrkN2++) {
1309 
1310  if(iTrkN2==iTrkP1 || iTrkN2==iTrkP2 || iTrkN2==iTrkN1) continue;
1311 
1312  //if(iTrkN2%1==0) AliDebug(1,Form(" 2nd loop on neg: track number %d of %d",iTrkN2,nSeleTrks));
1313 
1314  // get track from tracks array
1315  negtrack2 = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkN2);
1316 
1317  if(negtrack2->Charge()>0) continue;
1318 
1319  if(!TESTBIT(seleFlags[iTrkN2],kBitDispl)) continue;
1320 
1321  // Check single tracks cuts specific for 3 prongs
1322  if(!TESTBIT(seleFlags[iTrkN2],kBit3Prong)) continue;
1323  if(!TESTBIT(seleFlags[iTrkP1],kBit3Prong)) continue;
1324  if(!TESTBIT(seleFlags[iTrkN1],kBit3Prong)) continue;
1325 
1326  if(fMixEvent) {
1327  if(evtNumber[iTrkP1]==evtNumber[iTrkN2] ||
1328  evtNumber[iTrkN1]==evtNumber[iTrkN2] ||
1329  evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
1330  }
1331 
1332  if(isLikeSign2Prong) { // like-sign pair -> have to build only like-sign triplet
1333  if(!fLikeSign3prong) continue;
1334  if(postrack1->Charge()<0) { // ok: like-sign triplet (---)
1335  isLikeSign3Prong=kTRUE;
1336  } else { // not ok
1337  continue;
1338  }
1339  } else { // normal triplet (-+-)
1340  isLikeSign3Prong=kFALSE;
1341  }
1342 
1344  if(!TESTBIT(seleFlags[iTrkP1],kBitKaonCompat)) continue;
1345  }
1346  Bool_t okForLcTopKpi=kTRUE;
1347  Int_t pidLcStatus=3; // 3= OK as pKpi and Kpipi
1348  if(fUsePIDforLc>0){
1349  if(!TESTBIT(seleFlags[iTrkN1],kBitProtonCompat) &&
1350  !TESTBIT(seleFlags[iTrkN2],kBitProtonCompat) ){
1351  okForLcTopKpi=kFALSE;
1352  pidLcStatus=0;
1353  }
1354  if(okForLcTopKpi && fUsePIDforLc>1){
1355  okForLcTopKpi=kFALSE;
1356  pidLcStatus=0;
1357  if(TESTBIT(seleFlags[iTrkN1],kBitProtonCompat) &&
1358  TESTBIT(seleFlags[iTrkN2],kBitPionCompat) ){
1359  okForLcTopKpi=kTRUE;
1360  pidLcStatus+=1; // 1= OK as pKpi
1361  }
1362  if(TESTBIT(seleFlags[iTrkN2],kBitProtonCompat) &&
1363  TESTBIT(seleFlags[iTrkN1],kBitPionCompat) ){
1364  okForLcTopKpi=kTRUE;
1365  pidLcStatus+=2; // 2= OK as piKp
1366  }
1367  }
1368  }
1369  Bool_t okForDsToKKpi=kTRUE;
1370  if(fUseKaonPIDforDs){
1371  if(!TESTBIT(seleFlags[iTrkN1],kBitKaonCompat) &&
1372  !TESTBIT(seleFlags[iTrkN2],kBitKaonCompat) ) okForDsToKKpi=kFALSE;
1373  }
1374 
1375  // back to primary vertex
1376  // postrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
1377  // negtrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
1378  // negtrack2->PropagateToDCA(fV1,fBzkG,kVeryBig);
1379  SetParametersAtVertex(postrack1,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkP1));
1380  SetParametersAtVertex(negtrack1,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkN1));
1381  SetParametersAtVertex(negtrack2,(AliExternalTrackParam*)tracksAtVertex.UncheckedAt(iTrkN2));
1382  //printf("********** %d %d %d\n",postrack1->GetID(),negtrack1->GetID(),negtrack2->GetID());
1383 
1384  dcap1n2 = postrack1->GetDCA(negtrack2,fBzkG,xdummy,ydummy);
1385  if(dcap1n2>dcaMax) { negtrack2=0; continue; }
1386  dcan1n2 = negtrack1->GetDCA(negtrack2,fBzkG,xdummy,ydummy);
1387  if(dcan1n2>dcaMax) { negtrack2=0; continue; }
1388 
1389  threeTrackArray->AddAt(negtrack1,0);
1390  threeTrackArray->AddAt(postrack1,1);
1391  threeTrackArray->AddAt(negtrack2,2);
1392 
1393  // check invariant mass cuts for D+,Ds,Lc
1394  massCutOK=kTRUE;
1396  negtrack2->GetPxPyPz(momneg2);
1397  Double_t pxDau[3]={momneg1[0],mompos1[0],momneg2[0]};
1398  Double_t pyDau[3]={momneg1[1],mompos1[1],momneg2[1]};
1399  Double_t pzDau[3]={momneg1[2],mompos1[2],momneg2[2]};
1400  // massCutOK = SelectInvMassAndPt3prong(threeTrackArray);
1401  massCutOK = SelectInvMassAndPt3prong(pxDau,pyDau,pzDau,pidLcStatus);
1402  }
1403  if(!massCutOK) {
1404  threeTrackArray->Clear();
1405  negtrack2=0;
1406  continue;
1407  }
1408 
1409  // Vertexing
1410  twoTrackArray2->AddAt(postrack1,0);
1411  twoTrackArray2->AddAt(negtrack2,1);
1412 
1413  AliAODVertex *vertexp1n2 = ReconstructSecondaryVertex(twoTrackArray2,dispersion);
1414  if(!vertexp1n2) {
1415  twoTrackArray2->Clear();
1416  negtrack2=0;
1417  continue;
1418  }
1419 
1420  if(f3Prong) {
1421  AliAODVertex* secVert3PrAOD = ReconstructSecondaryVertex(threeTrackArray,dispersion);
1422  io3Prong = Make3Prong(threeTrackArray,event,secVert3PrAOD,dispersion,vertexp1n1,vertexp1n2,dcap1n1,dcap1n2,dcan1n2,okForLcTopKpi,okForDsToKKpi,ok3Prong);
1423  if(ok3Prong) {
1424  AliAODVertex *v3Prong = 0x0;
1425  if(!fMakeReducedRHF) v3Prong = new(verticesHFRef[iVerticesHF++])AliAODVertex(*secVert3PrAOD);
1426  if(!isLikeSign3Prong) {
1427  rd = new(aodCharm3ProngRef[i3Prong++])AliAODRecoDecayHF3Prong(*io3Prong);
1431  if(fMakeReducedRHF){
1432  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
1433  ((AliAODRecoDecayHF3Prong*)rd)->DeleteRecoD();
1434  }else{
1435  rd->SetSecondaryVtx(v3Prong);
1436  v3Prong->SetParent(rd);
1437  AddRefs(v3Prong,rd,event,threeTrackArray);
1438  }
1439  } else { // isLikeSign3Prong
1440  if(fLikeSign3prong){
1441  rd = new(aodLikeSign3ProngRef[iLikeSign3Prong++])AliAODRecoDecayHF3Prong(*io3Prong);
1445  if(fMakeReducedRHF){
1446  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
1447  ((AliAODRecoDecayHF3Prong*)rd)->DeleteRecoD();
1448  }else{
1449  rd->SetSecondaryVtx(v3Prong);
1450  v3Prong->SetParent(rd);
1451  AddRefs(v3Prong,rd,event,threeTrackArray);
1452  }
1453  }
1454 
1455  }
1456  }
1457  if(io3Prong) {delete io3Prong; io3Prong=NULL;}
1458  if(secVert3PrAOD) {delete secVert3PrAOD; secVert3PrAOD=NULL;}
1459  }
1460  threeTrackArray->Clear();
1461  negtrack2 = 0;
1462  delete vertexp1n2;
1463 
1464  } // end 2nd loop on negative tracks
1465 
1466  twoTrackArray2->Clear();
1467 
1468  negtrack1 = 0;
1469  delete vertexp1n1;
1470  } // end 1st loop on negative tracks
1471 
1472  postrack1 = 0;
1473  } // end 1st loop on positive tracks
1474 
1475 
1476  // AliDebug(1,Form(" Total HF vertices in event = %d;",
1477  // (Int_t)aodVerticesHFTClArr->GetEntriesFast()));
1478  if(fD0toKpi) {
1479  AliDebug(1,Form(" D0->Kpi in event = %d;",
1480  (Int_t)aodD0toKpiTClArr->GetEntriesFast()));
1481  }
1482  if(fJPSItoEle) {
1483  AliDebug(1,Form(" JPSI->ee in event = %d;",
1484  (Int_t)aodJPSItoEleTClArr->GetEntriesFast()));
1485  }
1486  if(f3Prong) {
1487  AliDebug(1,Form(" Charm->3Prong in event = %d;",
1488  (Int_t)aodCharm3ProngTClArr->GetEntriesFast()));
1489  }
1490  if(f4Prong) {
1491  AliDebug(1,Form(" Charm->4Prong in event = %d;\n",
1492  (Int_t)aodCharm4ProngTClArr->GetEntriesFast()));
1493  }
1494  if(fDstar) {
1495  AliDebug(1,Form(" D*->D0pi in event = %d;\n",
1496  (Int_t)aodDstarTClArr->GetEntriesFast()));
1497  }
1498  if(fCascades){
1499  AliDebug(1,Form(" cascades -> v0 + track in event = %d;\n",
1500  (Int_t)aodCascadesTClArr->GetEntriesFast()));
1501  }
1502  if(fLikeSign) {
1503  AliDebug(1,Form(" Like-sign 2Prong in event = %d;\n",
1504  (Int_t)aodLikeSign2ProngTClArr->GetEntriesFast()));
1505  }
1506  if(fLikeSign3prong && f3Prong) {
1507  AliDebug(1,Form(" Like-sign 3Prong in event = %d;\n",
1508  (Int_t)aodLikeSign3ProngTClArr->GetEntriesFast()));
1509  }
1510 
1511 
1512  twoTrackArray1->Delete(); delete twoTrackArray1;
1513  twoTrackArray2->Delete(); delete twoTrackArray2;
1514  twoTrackArrayCasc->Delete(); delete twoTrackArrayCasc;
1515  twoTrackArrayV0->Delete(); delete twoTrackArrayV0;
1516  threeTrackArray->Clear();
1517  threeTrackArray->Delete(); delete threeTrackArray;
1518  fourTrackArray->Delete(); delete fourTrackArray;
1519  delete [] seleFlags; seleFlags=NULL;
1520  if(evtNumber) {delete [] evtNumber; evtNumber=NULL;}
1521  tracksAtVertex.Delete();
1522 
1523  if(fInputAOD) {
1524  seleTrksArray.Delete();
1525  if(fAODMap) { delete [] fAODMap; fAODMap=NULL; }
1526  }
1527 
1528 
1529  //printf("Trks: total %d sele %d\n",fnTrksTotal,fnSeleTrksTotal);
1530 
1531  return;
1532 }
1533 //----------------------------------------------------------------------------
1535  const AliVEvent *event,
1536  const TObjArray *trkArray) const
1537 {
1540  //AliCodeTimerAuto("",0);
1541 
1542  if(fInputAOD) {
1543  AddDaughterRefs(v,event,trkArray);
1544  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
1545  }
1546 
1547  /*
1548  rd->SetListOfCutsRef((TList*)fListOfCuts);
1549  //fListOfCuts->Print();
1550  cout<<fListOfCuts<<endl;
1551  TList *l=(TList*)rd->GetListOfCuts();
1552  cout<<l<<endl;
1553  if(l) {l->Print(); }else{printf("error\n");}
1554  */
1555 
1556  return;
1557 }
1558 //---------------------------------------------------------------------------
1560  const AliVEvent *event,
1561  const TObjArray *trkArray) const
1562 {
1564  //AliCodeTimerAuto("",0);
1565 
1566  Int_t nDg = v->GetNDaughters();
1567  TObject *dg = 0;
1568  if(nDg) dg = v->GetDaughter(0);
1569 
1570  if(dg) return; // daughters already added
1571 
1572  Int_t nTrks = trkArray->GetEntriesFast();
1573 
1574  AliExternalTrackParam *track = 0;
1575  AliAODTrack *aodTrack = 0;
1576  Int_t id;
1577 
1578  for(Int_t i=0; i<nTrks; i++) {
1579  track = (AliExternalTrackParam*)trkArray->UncheckedAt(i);
1580  id = (Int_t)track->GetID();
1581  //printf("---> %d\n",id);
1582  if(id<0) continue; // this track is a AliAODRecoDecay
1583  aodTrack = dynamic_cast<AliAODTrack*>(event->GetTrack(fAODMap[id]));
1584  if(!aodTrack) AliFatal("Not a standard AOD");
1585  v->AddDaughter(aodTrack);
1586  }
1587 
1588  return;
1589 }
1590 //----------------------------------------------------------------------------
1592 {
1595  //
1596  //AliCodeTimerAuto("",0);
1597 
1598 
1599  TClonesArray *inputArray=(TClonesArray*)aod->GetList()->FindObject("VerticesHF");
1600  if(!inputArray) return;
1601 
1602  AliAODTrack *track = 0;
1603  AliAODVertex *vertex = 0;
1604 
1605  Bool_t needtofix=kFALSE;
1606  for(Int_t iv=0; iv<inputArray->GetEntriesFast(); iv++) {
1607  vertex = (AliAODVertex*)inputArray->UncheckedAt(iv);
1608  for(Int_t id=0; id<vertex->GetNDaughters(); id++) {
1609  track = (AliAODTrack*)vertex->GetDaughter(id);
1610  if(!track->GetStatus()) needtofix=kTRUE;
1611  }
1612  if(needtofix) break;
1613  }
1614 
1615  if(!needtofix) return;
1616 
1617 
1618  printf("Fixing references\n");
1619 
1620  fAODMapSize = 100000;
1621  fAODMap = new Int_t[fAODMapSize];
1622  memset(fAODMap,0,sizeof(Int_t)*fAODMapSize);
1623 
1624  for(Int_t i=0; i<aod->GetNumberOfTracks(); i++) {
1625  track = dynamic_cast<AliAODTrack*>(aod->GetTrack(i));
1626  if(!track) AliFatal("Not a standard AOD");
1627 
1628  // skip pure ITS SA tracks
1629  if(track->GetStatus()&AliESDtrack::kITSpureSA) continue;
1630 
1631  // skip tracks without ITS
1632  if(!(track->GetStatus()&AliESDtrack::kITSin)) continue;
1633 
1634  // TEMPORARY: check that the cov matrix is there
1635  Double_t covtest[21];
1636  if(!track->GetCovarianceXYZPxPyPz(covtest)) continue;
1637  //
1638 
1639  Int_t ind = (Int_t)track->GetID();
1640  if (ind>-1 && ind < fAODMapSize) fAODMap[ind] = i;
1641  }
1642 
1643 
1644  Int_t ids[4]={-1,-1,-1,-1};
1645  for(Int_t iv=0; iv<inputArray->GetEntriesFast(); iv++) {
1646  Bool_t cascade=kFALSE;
1647  vertex = (AliAODVertex*)inputArray->UncheckedAt(iv);
1648  Int_t id=0;
1649  Int_t nDgs = vertex->GetNDaughters();
1650  for(id=0; id<nDgs; id++) {
1651  track = (AliAODTrack*)vertex->GetDaughter(id);
1652  if(track->Charge()==0) {cascade=kTRUE; continue;} // cascade
1653  ids[id]=(Int_t)track->GetID();
1654  vertex->RemoveDaughter(track);
1655  }
1656  if(cascade) continue;
1657  for(id=0; id<nDgs; id++) {
1658  if (ids[id]>-1 && ids[id] < fAODMapSize) {
1659  track = dynamic_cast<AliAODTrack*>(aod->GetTrack(fAODMap[ids[id]]));
1660  if(!track) AliFatal("Not a standard AOD");
1661  vertex->AddDaughter(track);
1662  }
1663  }
1664 
1665  }
1666 
1667  return;
1668 }
1669 //----------------------------------------------------------------------------
1671  // method to retrieve daughters from trackID and reconstruct secondary vertex
1672  // save the TRefs to the candidate AliAODRecoDecayHF3Prong rd
1673  // and fill on-the-fly the data member of rd
1674  if(rd->GetIsFilled()!=0)return kTRUE;//if 0: reduced dAOD. skip if rd is already filled (1: standard dAOD, 2 already refilled)
1675  if(!fAODMap)MapAODtracks(event);//fill the AOD index map if it is not yet done
1676 
1677  TObjArray *threeTrackArray = new TObjArray(3);
1678 
1679  AliAODTrack *track1 =(AliAODTrack*)event->GetTrack(fAODMap[rd->GetProngID(0)]);//retrieve daughter from the trackID through the AOD index map
1680  if(!track1)return kFALSE;
1681  AliAODTrack *track2 =(AliAODTrack*)event->GetTrack(fAODMap[rd->GetProngID(1)]);
1682  if(!track2)return kFALSE;
1683  AliESDtrack *postrack1 = 0;
1684  AliESDtrack *negtrack1 = 0;
1685  postrack1 = new AliESDtrack(track1);
1686  negtrack1 = new AliESDtrack(track2);
1687 
1688  // DCA between the two tracks
1689  Double_t xdummy, ydummy;
1690  fBzkG = (Double_t)event->GetMagneticField();
1691  Double_t dca12 = postrack1->GetDCA(negtrack1,fBzkG,xdummy,ydummy);
1692 
1693  const AliVVertex *vprimary = event->GetPrimaryVertex();
1694  Double_t pos[3];
1695  Double_t cov[6];
1696  vprimary->GetXYZ(pos);
1697  vprimary->GetCovarianceMatrix(cov);
1698  fV1 = new AliESDVertex(pos,cov,100.,100,vprimary->GetName());
1699  fV1->GetCovMatrix(cov);
1700  if(!fVertexerTracks)fVertexerTracks=new AliVertexerTracks(fBzkG);
1701 
1702  AliAODTrack *track3 =(AliAODTrack*)event->GetTrack(fAODMap[rd->GetProngID(2)]);
1703  if(!track3)return kFALSE;
1704  AliESDtrack *esdt3 = new AliESDtrack(track3);
1705 
1706  Double_t dca2;
1707  Double_t dca3;
1708  threeTrackArray->AddAt(postrack1,0);
1709  threeTrackArray->AddAt(negtrack1,1);
1710  threeTrackArray->AddAt(esdt3,2);
1711  dca2 = esdt3->GetDCA(negtrack1,fBzkG,xdummy,ydummy);
1712  dca3 = esdt3->GetDCA(postrack1,fBzkG,xdummy,ydummy);
1713  Double_t dispersion;
1714 
1715  AliAODVertex* secVert3PrAOD = ReconstructSecondaryVertex(threeTrackArray, dispersion);
1716  if (!secVert3PrAOD) {
1717  threeTrackArray->Clear();
1718  threeTrackArray->Delete(); delete threeTrackArray;
1719  delete fV1; fV1=0;
1720  delete postrack1; postrack1=NULL;
1721  delete negtrack1; negtrack1=NULL;
1722  delete esdt3; esdt3=NULL;
1723  return kFALSE;
1724  }
1725 
1726  Double_t vtxRec=rd->GetDist12toPrim();
1727  Double_t vertexp2n1=rd->GetDist23toPrim();
1728  rd= Make3Prong(threeTrackArray, event, secVert3PrAOD,dispersion, vtxRec, vertexp2n1, dca12, dca2, dca3, rd);
1729  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
1730  rd->SetIsFilled(2);
1731  threeTrackArray->Clear();
1732  threeTrackArray->Delete(); delete threeTrackArray;
1733  delete fV1; fV1=0;
1734  delete postrack1; postrack1=NULL;
1735  delete negtrack1; negtrack1=NULL;
1736  delete esdt3; esdt3=NULL;
1737  return kTRUE;
1738 }
1739 //___________________________
1741  // method to retrieve daughters from trackID and reconstruct secondary vertex
1742  // save the TRefs to the candidate AliAODRecoDecayHF2Prong rd
1743  // and fill on-the-fly the data member of rd
1744  if(rd->GetIsFilled()!=0)return kTRUE;//if 0: reduced dAOD. skip if rd is already filled (1:standard dAOD, 2 already refilled)
1745  if(!fAODMap)MapAODtracks(event);//fill the AOD index map if it is not yet done
1746 
1747  Double_t dispersion;
1748  TObjArray *twoTrackArray1 = new TObjArray(2);
1749 
1750  AliAODTrack *track1 =(AliAODTrack*)event->GetTrack(fAODMap[rd->GetProngID(0)]);//retrieve daughter from the trackID through the AOD index map
1751  if(!track1)return kFALSE;
1752  AliAODTrack *track2 =(AliAODTrack*)event->GetTrack(fAODMap[rd->GetProngID(1)]);
1753  if(!track2)return kFALSE;
1754 
1755  AliESDtrack *esdt1 = 0;
1756  AliESDtrack *esdt2 = 0;
1757  esdt1 = new AliESDtrack(track1);
1758  esdt2 = new AliESDtrack(track2);
1759 
1760  twoTrackArray1->AddAt(esdt1,0);
1761  twoTrackArray1->AddAt(esdt2,1);
1762  // DCA between the two tracks
1763  Double_t xdummy, ydummy;
1764  fBzkG = (Double_t)event->GetMagneticField();
1765  Double_t dca12 = esdt1->GetDCA(esdt2,fBzkG,xdummy,ydummy);
1766  const AliVVertex *vprimary = event->GetPrimaryVertex();
1767  Double_t pos[3];
1768  Double_t cov[6];
1769  vprimary->GetXYZ(pos);
1770  vprimary->GetCovarianceMatrix(cov);
1771  fV1 = new AliESDVertex(pos,cov,100.,100,vprimary->GetName());
1772  fV1->GetCovMatrix(cov);
1773  if(!fVertexerTracks)fVertexerTracks=new AliVertexerTracks(fBzkG);
1774 
1775 
1776  AliAODVertex *vtxRec = ReconstructSecondaryVertex(twoTrackArray1, dispersion);
1777  if(!vtxRec) {
1778  twoTrackArray1->Clear();
1779  twoTrackArray1->Delete(); delete twoTrackArray1;
1780  delete fV1; fV1=0;
1781  delete esdt1; esdt1=NULL;
1782  delete esdt2; esdt2=NULL;
1783  return kFALSE; }
1784  Bool_t okD0=kFALSE;
1785  Bool_t okJPSI=kFALSE;
1786  Bool_t okD0FromDstar=kFALSE;
1787  Bool_t refill =kTRUE;
1788  rd= Make2Prong(twoTrackArray1, event, vtxRec, dca12, okD0, okJPSI, okD0FromDstar,refill,rd);
1789  rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
1790  rd->SetIsFilled(2);
1791  delete fV1; fV1=0;
1792  twoTrackArray1->Clear();
1793  twoTrackArray1->Delete(); delete twoTrackArray1;
1794  delete esdt1; esdt1=NULL;
1795  delete esdt2; esdt2=NULL;
1796  return kTRUE;
1797 }
1798 //----------------------------------------------------------------------------
1799 Bool_t AliAnalysisVertexingHF::FillRecoCasc(AliVEvent *event,AliAODRecoCascadeHF *rCasc, Bool_t DStar){
1800  // method to retrieve daughters from trackID
1801  // and fill on-the-fly the data member of rCasc and their AliAODRecoDecayHF2Prong daughters
1802  if(rCasc->GetIsFilled()!=0) return kTRUE;//if 0: reduced dAOD. skip if rd is already filled (1: standard dAOD, 2: already refilled)
1803  if(!fAODMap)MapAODtracks(event);//fill the AOD index map if it is not yet done
1804  TObjArray *twoTrackArrayCasc = new TObjArray(2);
1805 
1806  AliAODTrack *trackB =(AliAODTrack*)event->GetTrack(fAODMap[rCasc->GetProngID(0)]);//retrieve bachelor from the trackID through the AOD index map
1807  if(!trackB)return kFALSE;
1808 
1809  AliNeutralTrackParam *trackV0=NULL;
1810  AliAODv0 *v0 =NULL;
1811  AliAODRecoDecayHF2Prong *trackD0=NULL;
1812 
1813  if(DStar){
1814  TClonesArray *inputArrayD0=(TClonesArray*)event->GetList()->FindObject("D0toKpi");
1815  if(!inputArrayD0) return kFALSE;
1816  trackD0=(AliAODRecoDecayHF2Prong*)inputArrayD0->At(rCasc->GetProngID(1));
1817  if(!trackD0)return kFALSE;
1818  FillRecoCand(event,trackD0);//fill missing info of the corresponding D0 daughter
1819 
1820  trackV0 = new AliNeutralTrackParam(trackD0);
1821 
1822  }else{//is a V0 candidate
1823  v0 = ((AliAODEvent*)event)->GetV0(rCasc->GetProngID(1));
1824  if(!v0) return kFALSE;
1825  // Define the V0 (neutral) track
1826  const AliVTrack *trackVV0 = dynamic_cast<const AliVTrack*>(v0);
1827  if(trackVV0) trackV0 = new AliNeutralTrackParam(trackVV0);
1828  }
1829 
1830  AliESDtrack *esdB = new AliESDtrack(trackB);
1831 
1832  twoTrackArrayCasc->AddAt(esdB,0);
1833  twoTrackArrayCasc->AddAt(trackV0,1);
1834 
1835  fBzkG = (Double_t)event->GetMagneticField();
1836  const AliVVertex *vprimary = event->GetPrimaryVertex();
1837 
1838  Double_t pos[3];
1839  Double_t cov[6];
1840  vprimary->GetXYZ(pos);
1841  vprimary->GetCovarianceMatrix(cov);
1842  fV1 = new AliESDVertex(pos,cov,100.,100,vprimary->GetName());
1843  fV1->GetCovMatrix(cov);
1844 
1845  Double_t dca = 0.;
1846  AliAODVertex *vtxCasc = 0x0;
1847  Double_t chi2perNDF = fV1->GetChi2toNDF();
1848  vtxCasc = new AliAODVertex(pos,cov,chi2perNDF,0x0,-1,AliAODVertex::kUndef,2);
1849  if(!vtxCasc) {
1850  twoTrackArrayCasc->Clear();
1851  twoTrackArrayCasc->Delete(); delete twoTrackArrayCasc;
1852  delete fV1; fV1=0;
1853  delete esdB; esdB=NULL;
1854  delete vtxCasc;vtxCasc=NULL;
1855  delete trackB; trackB=NULL;
1856  delete trackV0; trackV0=NULL;
1857  if(!DStar){
1858  v0=NULL;
1859  }
1860  return kFALSE;
1861  }
1862  vtxCasc->SetParent(rCasc);
1863  rCasc->SetSecondaryVtx(vtxCasc);
1864  AddDaughterRefs(vtxCasc,(AliAODEvent*)event,twoTrackArrayCasc);
1865  if(DStar)vtxCasc->AddDaughter(trackD0);
1866  else vtxCasc->AddDaughter(v0);
1867  rCasc->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
1868 
1869  Bool_t refill =kTRUE;
1870  Bool_t dummy1, dummy2, dummy3;
1871 
1872  Double_t px[2],py[2],pz[2],d0[2],d0err[2];
1873  // propagate tracks to secondary vertex, to compute inv. mass
1874  esdB->PropagateToDCA(vtxCasc,fBzkG,kVeryBig);
1875  trackV0->PropagateToDCA(vtxCasc,fBzkG,kVeryBig);
1876  Double_t momentum[3];
1877  esdB->GetPxPyPz(momentum);
1878  px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2];
1879  trackV0->GetPxPyPz(momentum);
1880  px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2];
1881 
1882  AliAODVertex *primVertexAOD = PrimaryVertex(twoTrackArrayCasc,event);
1883  if(!primVertexAOD){
1884  delete fV1; fV1=0;
1885  delete vtxCasc; vtxCasc=NULL;
1886  twoTrackArrayCasc->Clear();
1887  twoTrackArrayCasc->Delete(); delete twoTrackArrayCasc;
1888  delete esdB; esdB=NULL;
1889  delete trackV0; trackV0=NULL;
1890  if(!DStar)v0=NULL;
1891  return kFALSE;
1892  }
1893  Double_t d0z0[2],covd0z0[3];
1894  esdB->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
1895  d0[0] = d0z0[0];
1896  d0err[0] = TMath::Sqrt(covd0z0[0]);
1897  trackV0->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
1898  d0[1] = d0z0[0];
1899  d0err[1] = TMath::Sqrt(covd0z0[0]);
1900  rCasc->SetPxPyPzProngs(2,px,py,pz);
1901  rCasc->SetDCA(dca);
1902  rCasc->Setd0Prongs(2,d0);
1903  rCasc->Setd0errProngs(2,d0err);
1904  rCasc->SetCharge(esdB->Charge());
1905  // get PID info from ESD
1906  Double_t esdpid0[5]={0.,0.,0.,0.,0.};
1907  if(esdB->GetStatus()&AliESDtrack::kESDpid) esdB->GetESDpid(esdpid0);
1908  Double_t esdpid1[5]={0.,0.,0.,0.,0.};
1909  Double_t esdpid[10];
1910  for(Int_t i=0;i<5;i++) {
1911  esdpid[i] = esdpid0[i];
1912  esdpid[5+i] = esdpid1[i];
1913  }
1914  rCasc->SetPID(2,esdpid);
1915  rCasc->SetIsFilled(2);
1916 
1917 
1918  delete fV1; fV1=0;
1919  if(primVertexAOD) {delete primVertexAOD; primVertexAOD=NULL;}
1920  twoTrackArrayCasc->Clear();
1921  twoTrackArrayCasc->Delete(); delete twoTrackArrayCasc;
1922  delete esdB; esdB=NULL;
1923  delete trackV0; trackV0=NULL;
1924  if(!DStar)v0=NULL;
1925  return kTRUE;
1926 }
1927 //---------------------------------------------------------------------------
1929  TObjArray *twoTrackArray,AliVEvent *event,
1930  AliAODVertex *secVert,
1931  AliAODRecoDecayHF2Prong *rd2Prong,
1932  Double_t dca,
1933  Bool_t &okDstar)
1934 {
1937  //AliCodeTimerAuto("",0);
1938  UInt_t ntref=TProcessID::GetObjectCount();
1939  if(ntref>16776216){// This number is 2^24-1000. The maximum number of TRef for a given TProcesssID is 2^24=16777216.
1940  AliFatal(Form("Number of TRef created (=%d), close to limit (16777216)",ntref));
1941  }
1942 
1943  okDstar = kFALSE;
1944 
1945  Bool_t dummy1,dummy2,dummy3;
1946 
1947  // We use Make2Prong to construct the AliAODRecoCascadeHF
1948  // (which inherits from AliAODRecoDecayHF2Prong)
1949  AliAODRecoCascadeHF *theCascade =
1950  (AliAODRecoCascadeHF*)Make2Prong(twoTrackArray,event,secVert,dca,
1951  dummy1,dummy2,dummy3);
1952  if(!theCascade) return 0x0;
1953 
1954  // charge
1955  AliESDtrack *trackPi = (AliESDtrack*)twoTrackArray->UncheckedAt(0);
1956  theCascade->SetCharge(trackPi->Charge());
1957 
1958  //--- selection cuts
1959  //
1960  AliAODRecoCascadeHF *tmpCascade = new AliAODRecoCascadeHF(*theCascade);
1961  if(fInputAOD){
1962  Int_t idSoftPi=(Int_t)trackPi->GetID();
1963  if (idSoftPi > -1 && idSoftPi < fAODMapSize) {
1964  AliAODTrack* trackPiAOD=dynamic_cast<AliAODTrack*>(event->GetTrack(fAODMap[idSoftPi]));
1965  if(!trackPiAOD) AliFatal("Not a standard AOD");
1966  tmpCascade->GetSecondaryVtx()->AddDaughter(trackPiAOD);
1967  }
1968  }else{
1969  tmpCascade->GetSecondaryVtx()->AddDaughter(trackPi);
1970  }
1971  tmpCascade->GetSecondaryVtx()->AddDaughter(rd2Prong);
1972 
1973  AliAODVertex *primVertexAOD=0;
1975  // take event primary vertex
1976  primVertexAOD = PrimaryVertex();
1977  tmpCascade->SetOwnPrimaryVtx(primVertexAOD);
1978  rd2Prong->SetOwnPrimaryVtx(primVertexAOD);
1979  }
1980  // select D*->D0pi
1981  if(fDstar) {
1982  okDstar = (Bool_t)fCutsDStartoKpipi->IsSelected(tmpCascade,AliRDHFCuts::kCandidate);
1983  if(okDstar) theCascade->SetSelectionBit(AliRDHFCuts::kDstarCuts);
1984  }
1985  tmpCascade->GetSecondaryVtx()->RemoveDaughters();
1986  tmpCascade->UnsetOwnPrimaryVtx();
1987  delete tmpCascade; tmpCascade=NULL;
1989  rd2Prong->UnsetOwnPrimaryVtx();
1990  }
1991  if(primVertexAOD) {delete primVertexAOD; primVertexAOD=NULL;}
1992  //---
1993 
1994 
1995  return theCascade;
1996 }
1997 
1998 
1999 //----------------------------------------------------------------------------
2001  TObjArray *twoTrackArray,AliVEvent *event,
2002  AliAODVertex *secVert,
2003  AliAODv0 *v0,
2004  Double_t dca,
2005  Bool_t &okCascades)
2006 {
2009  //AliCodeTimerAuto("",0);
2010  UInt_t ntref=TProcessID::GetObjectCount();
2011  if(ntref>16776216){// This number is 2^24-1000. The maximum number of TRef for a given TProcesssID is 2^24=16777216.
2012  AliFatal(Form("Number of TRef created (=%d), close to limit (16777216)",ntref));
2013  }
2014 
2015  // AliDebug(2,Form(" building the cascade"));
2016  okCascades= kFALSE;
2017  Bool_t dummy1,dummy2,dummy3;
2018 
2019  // We use Make2Prong to construct the AliAODRecoCascadeHF
2020  // (which inherits from AliAODRecoDecayHF2Prong)
2021  AliAODRecoCascadeHF *theCascade =
2022  (AliAODRecoCascadeHF*)Make2Prong(twoTrackArray,event,secVert,dca,
2023  dummy1,dummy2,dummy3);
2024  if(!theCascade) return 0x0;
2025 
2026  // bachelor track and charge
2027  AliESDtrack *trackBachelor = (AliESDtrack*)twoTrackArray->UncheckedAt(0);
2028  theCascade->SetCharge(trackBachelor->Charge());
2029 
2030  //--- selection cuts
2031  //
2032 
2033  AliAODRecoCascadeHF *tmpCascade = new AliAODRecoCascadeHF(*theCascade);
2034  if(fInputAOD){
2035  Int_t idBachelor=(Int_t)trackBachelor->GetID();
2036  if (idBachelor > -1 && idBachelor < fAODMapSize) {
2037  AliAODTrack* trackBachelorAOD=dynamic_cast<AliAODTrack*>(event->GetTrack(fAODMap[idBachelor]));
2038  if(!trackBachelorAOD) AliFatal("Not a standard AOD");
2039  tmpCascade->GetSecondaryVtx()->AddDaughter(trackBachelorAOD);
2040  }
2041  }else{
2042  tmpCascade->GetSecondaryVtx()->AddDaughter(trackBachelor);
2043  }
2044  tmpCascade->GetSecondaryVtx()->AddDaughter(v0);
2045 
2046  AliAODVertex *primVertexAOD=0;
2048  // take event primary vertex
2049  primVertexAOD = PrimaryVertex();
2050  if(!primVertexAOD) primVertexAOD = (AliAODVertex*)event->GetPrimaryVertex();
2051  tmpCascade->SetOwnPrimaryVtx(primVertexAOD);
2052  }
2053 
2054  // select Cascades
2055  if(fCascades && fInputAOD){
2056  okCascades = (Bool_t)fCutsLctoV0->IsSelected(tmpCascade,AliRDHFCuts::kCandidate);
2057  }
2058  else {
2059  //AliDebug(2,Form("The cascade is contructed from ESDs, no cuts are applied"));
2060  okCascades=kTRUE;
2061  }// no cuts implemented from ESDs
2062  tmpCascade->GetSecondaryVtx()->RemoveDaughters();
2063  tmpCascade->UnsetOwnPrimaryVtx();
2064  delete tmpCascade; tmpCascade=NULL;
2065  if(primVertexAOD) {delete primVertexAOD; primVertexAOD=NULL;}
2066  //---
2067 
2068  return theCascade;
2069 }
2070 
2071 //-----------------------------------------------------------------------------
2073  TObjArray *twoTrackArray,AliVEvent *event,
2074  AliAODVertex *secVert,Double_t dca,
2075  Bool_t &okD0,Bool_t &okJPSI,
2076  Bool_t &okD0fromDstar, Bool_t refill, AliAODRecoDecayHF2Prong *rd)
2077 {
2080  // G.E.Bruno (J/psi), A.Dainese (D0->Kpi)
2081  // AliCodeTimerAuto("",0);
2082  okD0=kFALSE; okJPSI=kFALSE; okD0fromDstar=kFALSE;
2083 
2084  UInt_t ntref=TProcessID::GetObjectCount();
2085  if(ntref>16776216){// This number is 2^24-1000. The maximum number of TRef for a given TProcesssID is 2^24=16777216.
2086  AliFatal(Form("Number of TRef created (=%d), close to limit (16777216)",ntref));
2087  }
2088 
2089  Double_t px[2],py[2],pz[2],d0[2],d0err[2];
2090  AliESDtrack *postrack = (AliESDtrack*)twoTrackArray->UncheckedAt(0);
2091  AliESDtrack *negtrack = (AliESDtrack*)twoTrackArray->UncheckedAt(1);
2092 
2093  // propagate tracks to secondary vertex, to compute inv. mass
2094  postrack->PropagateToDCA(secVert,fBzkG,kVeryBig);
2095  negtrack->PropagateToDCA(secVert,fBzkG,kVeryBig);
2096 
2097  Double_t momentum[3];
2098  postrack->GetPxPyPz(momentum);
2099  px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2];
2100  negtrack->GetPxPyPz(momentum);
2101  px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2];
2102 
2103  if(!refill){//skip if it is called in refill step because already checked
2104  // invariant mass cut (try to improve coding here..)
2105  Bool_t okMassCut=kFALSE;
2106  if(!okMassCut && fD0toKpi) if(SelectInvMassAndPtD0Kpi(px,py,pz)) okMassCut=kTRUE;
2107  if(!okMassCut && fJPSItoEle) if(SelectInvMassAndPtJpsiee(px,py,pz)) okMassCut=kTRUE;
2108  if(!okMassCut && fDstar) if(SelectInvMassAndPtDstarD0pi(px,py,pz)) okMassCut=kTRUE;
2109  if(!okMassCut && fCascades) if(SelectInvMassAndPtCascade(px,py,pz)) okMassCut=kTRUE;
2110  if(!okMassCut) {
2111  //AliDebug(2," candidate didn't pass mass cut");
2112  return 0x0;
2113  }
2114  }
2115  // primary vertex to be used by this candidate
2116  AliAODVertex *primVertexAOD = PrimaryVertex(twoTrackArray,event);
2117  if(!primVertexAOD) return 0x0;
2118 
2119  Double_t d0z0[2],covd0z0[3];
2120  postrack->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2121  d0[0] = d0z0[0];
2122  d0err[0] = TMath::Sqrt(covd0z0[0]);
2123  negtrack->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2124  d0[1] = d0z0[0];
2125  d0err[1] = TMath::Sqrt(covd0z0[0]);
2126  AliAODRecoDecayHF2Prong *the2Prong;
2127  // create the object AliAODRecoDecayHF2Prong
2128  if(!refill){
2129  the2Prong = new AliAODRecoDecayHF2Prong(secVert,px,py,pz,d0,d0err,dca);
2130  the2Prong->SetOwnPrimaryVtx(primVertexAOD);
2131  UShort_t id[2]={(UShort_t)postrack->GetID(),(UShort_t)negtrack->GetID()};
2132  the2Prong->SetProngIDs(2,id);
2133  if(postrack->Charge()!=0 && negtrack->Charge()!=0) { // don't apply these cuts if it's a Dstar
2134  // Add daughter references already here
2135  if(fInputAOD) AddDaughterRefs(secVert,(AliAODEvent*)event,twoTrackArray);
2136 
2137  // select D0->Kpi
2138  if(fD0toKpi) {
2139  okD0 = (Bool_t)fCutsD0toKpi->IsSelected(the2Prong,AliRDHFCuts::kCandidate,(AliAODEvent*)event);
2140  if(okD0) the2Prong->SetSelectionBit(AliRDHFCuts::kD0toKpiCuts);
2141  }
2142  //if(fDebug && fD0toKpi) printf(" %d\n",(Int_t)okD0);
2143  // select J/psi from B
2144  if(fJPSItoEle) {
2145  okJPSI = (Bool_t)fCutsJpsitoee->IsSelected(the2Prong,AliRDHFCuts::kCandidate);
2146  }
2147  //if(fDebug && fJPSItoEle) printf(" %d\n",(Int_t)okJPSI);
2148  // select D0->Kpi from Dstar
2149  if(fDstar) {
2150  okD0fromDstar = (Bool_t)fCutsDStartoKpipi->IsD0FromDStarSelected(the2Prong->Pt(),the2Prong,AliRDHFCuts::kCandidate);
2151  if(okD0fromDstar) the2Prong->SetSelectionBit(AliRDHFCuts::kD0fromDstarCuts);
2152  }
2153  //if(fDebug && fDstar) printf(" %d\n",(Int_t)okD0fromDstar);
2154  }
2155  }else{
2156  the2Prong =rd;
2157  the2Prong->SetSecondaryVtx(secVert);
2158  secVert->SetParent(the2Prong);
2159  AddDaughterRefs(secVert,(AliAODEvent*)event,twoTrackArray);
2160  the2Prong->SetOwnPrimaryVtx(primVertexAOD);
2161  the2Prong->SetPxPyPzProngs(2,px,py,pz);
2162  the2Prong->SetDCA(dca);
2163  the2Prong->Setd0Prongs(2,d0);
2164  the2Prong->Setd0errProngs(2,d0err);
2165  the2Prong->SetCharge(0);
2166  }
2167  delete primVertexAOD; primVertexAOD=NULL;
2168 
2169  // remove the primary vertex (was used only for selection)
2171  the2Prong->UnsetOwnPrimaryVtx();
2172  }
2173 
2174  // get PID info from ESD
2175  Double_t esdpid0[5]={0.,0.,0.,0.,0.};
2176  if(postrack->GetStatus()&AliESDtrack::kESDpid) postrack->GetESDpid(esdpid0);
2177  Double_t esdpid1[5]={0.,0.,0.,0.,0.};
2178  if(negtrack->GetStatus()&AliESDtrack::kESDpid) negtrack->GetESDpid(esdpid1);
2179  Double_t esdpid[10];
2180  for(Int_t i=0;i<5;i++) {
2181  esdpid[i] = esdpid0[i];
2182  esdpid[5+i] = esdpid1[i];
2183  }
2184  the2Prong->SetPID(2,esdpid);
2185  return the2Prong;
2186 }
2187 //----------------------------------------------------------------------------
2189  TObjArray *threeTrackArray,AliVEvent *event,
2190  AliAODVertex *secVert,Double_t dispersion,
2191  const AliAODVertex *vertexp1n1,const AliAODVertex *vertexp2n1,
2192  Double_t dcap1n1,Double_t dcap2n1,Double_t dcap1p2,
2193  Bool_t useForLc, Bool_t useForDs, Bool_t &ok3Prong)
2194 {
2197  // E.Bruna, F.Prino
2198  // AliCodeTimerAuto("",0);
2199 
2200 
2201  UInt_t ntref=TProcessID::GetObjectCount();
2202  if(ntref>16776216){// This number is 2^24-1000. The maximum number of TRef for a given TProcesssID is 2^24=16777216.
2203  AliFatal(Form("Number of TRef created (=%d), close to limit (16777216)",ntref));
2204  }
2205 
2206  ok3Prong=kFALSE;
2207  if(!secVert || !vertexp1n1 || !vertexp2n1) return 0x0;
2208 
2209  Double_t px[3],py[3],pz[3],d0[3],d0err[3];
2210  Double_t momentum[3];
2211 
2212 
2213  AliESDtrack *postrack1 = (AliESDtrack*)threeTrackArray->UncheckedAt(0);
2214  AliESDtrack *negtrack = (AliESDtrack*)threeTrackArray->UncheckedAt(1);
2215  AliESDtrack *postrack2 = (AliESDtrack*)threeTrackArray->UncheckedAt(2);
2216 
2217  postrack1->PropagateToDCA(secVert,fBzkG,kVeryBig);
2218  negtrack->PropagateToDCA(secVert,fBzkG,kVeryBig);
2219  postrack2->PropagateToDCA(secVert,fBzkG,kVeryBig);
2220  postrack1->GetPxPyPz(momentum);
2221  px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2];
2222  negtrack->GetPxPyPz(momentum);
2223  px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2];
2224  postrack2->GetPxPyPz(momentum);
2225  px[2] = momentum[0]; py[2] = momentum[1]; pz[2] = momentum[2];
2226 
2227  // invariant mass cut for D+, Ds, Lc
2228  Bool_t okMassCut=kFALSE;
2229  if(fMassCutBeforeVertexing) okMassCut=kTRUE; // mass cut already done and passed
2230  if(!okMassCut && f3Prong) if(SelectInvMassAndPt3prong(px,py,pz)) okMassCut=kTRUE;
2231  if(!okMassCut) {
2232  //AliDebug(2," candidate didn't pass mass cut");
2233  return 0x0;
2234  }
2235 
2236  // primary vertex to be used by this candidate
2237  AliAODVertex *primVertexAOD = PrimaryVertex(threeTrackArray,event);
2238  if(!primVertexAOD) return 0x0;
2239 
2240  Double_t d0z0[2],covd0z0[3];
2241  postrack1->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2242  d0[0]=d0z0[0];
2243  d0err[0] = TMath::Sqrt(covd0z0[0]);
2244  negtrack->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2245  d0[1]=d0z0[0];
2246  d0err[1] = TMath::Sqrt(covd0z0[0]);
2247  postrack2->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2248  d0[2]=d0z0[0];
2249  d0err[2] = TMath::Sqrt(covd0z0[0]);
2250 
2251 
2252  // create the object AliAODRecoDecayHF3Prong
2253  Double_t pos[3]; primVertexAOD->GetXYZ(pos);
2254  Double_t dca[3]={dcap1n1,dcap2n1,dcap1p2};
2255  Double_t dist12=TMath::Sqrt((vertexp1n1->GetX()-pos[0])*(vertexp1n1->GetX()-pos[0])+(vertexp1n1->GetY()-pos[1])*(vertexp1n1->GetY()-pos[1])+(vertexp1n1->GetZ()-pos[2])*(vertexp1n1->GetZ()-pos[2]));
2256  Double_t dist23=TMath::Sqrt((vertexp2n1->GetX()-pos[0])*(vertexp2n1->GetX()-pos[0])+(vertexp2n1->GetY()-pos[1])*(vertexp2n1->GetY()-pos[1])+(vertexp2n1->GetZ()-pos[2])*(vertexp2n1->GetZ()-pos[2]));
2257  Short_t charge=(Short_t)(postrack1->Charge()+postrack2->Charge()+negtrack->Charge());
2258  AliAODRecoDecayHF3Prong *the3Prong = new AliAODRecoDecayHF3Prong(secVert,px,py,pz,d0,d0err,dca,dispersion,dist12,dist23,charge);
2259  the3Prong->SetOwnPrimaryVtx(primVertexAOD);
2260  UShort_t id[3]={(UShort_t)postrack1->GetID(),(UShort_t)negtrack->GetID(),(UShort_t)postrack2->GetID()};
2261  the3Prong->SetProngIDs(3,id);
2262 
2263  delete primVertexAOD; primVertexAOD=NULL;
2264 
2265  // Add daughter references already here
2266  if(fInputAOD) AddDaughterRefs(secVert,(AliAODEvent*)event,threeTrackArray);
2267 
2268  // select D+->Kpipi, Ds->KKpi, Lc->pKpi
2269  if(f3Prong) {
2270  ok3Prong = kFALSE;
2271 
2272  if(fOKInvMassDplus && fCutsDplustoKpipi->IsSelected(the3Prong,AliRDHFCuts::kCandidate,(AliAODEvent*)event)) {
2273  ok3Prong = kTRUE;
2275  }
2276  if(useForDs && fOKInvMassDs){
2277  if(fCutsDstoKKpi->IsSelected(the3Prong,AliRDHFCuts::kCandidate,(AliAODEvent*)event)) {
2278  ok3Prong = kTRUE;
2280  }
2281  }
2282  if(useForLc && fOKInvMassLc){
2283  if(fCutsLctopKpi->IsSelected(the3Prong,AliRDHFCuts::kCandidate,(AliAODEvent*)event)) {
2284  ok3Prong = kTRUE;
2286  }
2287  }
2288  }
2289  //if(fDebug) printf("ok3Prong: %d\n",(Int_t)ok3Prong);
2290 
2292  the3Prong->UnsetOwnPrimaryVtx();
2293  }
2294 
2295  // get PID info from ESD
2296  Double_t esdpid0[5]={0.,0.,0.,0.,0.};
2297  if(postrack1->GetStatus()&AliESDtrack::kESDpid) postrack1->GetESDpid(esdpid0);
2298  Double_t esdpid1[5]={0.,0.,0.,0.,0.};
2299  if(negtrack->GetStatus()&AliESDtrack::kESDpid) negtrack->GetESDpid(esdpid1);
2300  Double_t esdpid2[5]={0.,0.,0.,0.,0.};
2301  if(postrack2->GetStatus()&AliESDtrack::kESDpid) postrack2->GetESDpid(esdpid2);
2302 
2303  Double_t esdpid[15];
2304  for(Int_t i=0;i<5;i++) {
2305  esdpid[i] = esdpid0[i];
2306  esdpid[5+i] = esdpid1[i];
2307  esdpid[10+i] = esdpid2[i];
2308  }
2309  the3Prong->SetPID(3,esdpid);
2310 
2311  return the3Prong;
2312 }
2313 //----------------------------------------------------------------------------
2315  TObjArray *threeTrackArray,AliVEvent *event,
2316  AliAODVertex *secVert,Double_t dispersion,
2317  Double32_t dist12, Double32_t dist23,
2318  Double_t dcap1n1,Double_t dcap2n1,Double_t dcap1p2,
2320 {
2321  // Fill on-the-fly the 3prong data member missing info
2322  // set TRef of vertex and daughters
2323  // do not recalculate the two-track secondary vertex dist12 and dist23
2324  // because they are stored in the AliAODRecoDecayHF3Prong candidate
2325  // reconstructed in the FindCandidate step
2326  // do not check if it is a Lambdac, Ds or Dplus because it is already check in the FindCandidate step
2327  // and the info is stored
2328  // AliCodeTimerAuto("",0);
2329 
2330  UInt_t ntref=TProcessID::GetObjectCount();
2331  if(ntref>16776216){// This number is 2^24-1000. The maximum number of TRef for a given TProcesssID is 2^24=16777216.
2332  AliFatal(Form("Number of TRef created (=%d), close to limit (16777216)",ntref));
2333  }
2334 
2335  Double_t px[3],py[3],pz[3],d0[3],d0err[3];
2336  Double_t momentum[3];
2337 
2338 
2339  AliESDtrack *postrack1 = (AliESDtrack*)threeTrackArray->UncheckedAt(0);
2340  AliESDtrack *negtrack = (AliESDtrack*)threeTrackArray->UncheckedAt(1);
2341  AliESDtrack *postrack2 = (AliESDtrack*)threeTrackArray->UncheckedAt(2);
2342 
2343  postrack1->PropagateToDCA(secVert,fBzkG,kVeryBig);
2344  negtrack->PropagateToDCA(secVert,fBzkG,kVeryBig);
2345  postrack2->PropagateToDCA(secVert,fBzkG,kVeryBig);
2346  postrack1->GetPxPyPz(momentum);
2347  px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2];
2348  negtrack->GetPxPyPz(momentum);
2349  px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2];
2350  postrack2->GetPxPyPz(momentum);
2351  px[2] = momentum[0]; py[2] = momentum[1]; pz[2] = momentum[2];
2352  // primary vertex to be used by this candidate
2353  AliAODVertex *primVertexAOD = PrimaryVertex(threeTrackArray,event);
2354  if(!primVertexAOD) return 0x0;
2355  Double_t d0z0[2],covd0z0[3];
2356  postrack1->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2357  d0[0]=d0z0[0];
2358  d0err[0] = TMath::Sqrt(covd0z0[0]);
2359  negtrack->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2360  d0[1]=d0z0[0];
2361  d0err[1] = TMath::Sqrt(covd0z0[0]);
2362  postrack2->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2363  d0[2]=d0z0[0];
2364  d0err[2] = TMath::Sqrt(covd0z0[0]);
2365 
2366  Double_t pos[3]; primVertexAOD->GetXYZ(pos);
2367  Double_t dca[3]={dcap1n1,dcap2n1,dcap1p2};
2368  Short_t charge=(Short_t)(postrack1->Charge()+postrack2->Charge()+negtrack->Charge());
2369 
2370  rd->SetSecondaryVtx(secVert);
2371  secVert->SetParent(rd);
2372  AddDaughterRefs(secVert,(AliAODEvent*)event,threeTrackArray);
2373 
2374  rd->SetPxPyPzProngs(3,px,py,pz);
2375  rd->SetDCAs(3,dca);
2376  rd->Setd0Prongs(3,d0);
2377  rd->Setd0errProngs(3,d0err);
2378  rd->SetCharge(charge);
2379  rd->SetOwnPrimaryVtx(primVertexAOD);
2380  rd->SetSigmaVert(dispersion);
2381  delete primVertexAOD; primVertexAOD=NULL;
2382 
2384  rd->UnsetOwnPrimaryVtx();
2385  }
2386 
2387  // get PID info from ESD
2388  Double_t esdpid0[5]={0.,0.,0.,0.,0.};
2389  if(postrack1->GetStatus()&AliESDtrack::kESDpid) postrack1->GetESDpid(esdpid0);
2390  Double_t esdpid1[5]={0.,0.,0.,0.,0.};
2391  if(negtrack->GetStatus()&AliESDtrack::kESDpid) negtrack->GetESDpid(esdpid1);
2392  Double_t esdpid2[5]={0.,0.,0.,0.,0.};
2393  if(postrack2->GetStatus()&AliESDtrack::kESDpid) postrack2->GetESDpid(esdpid2);
2394 
2395  Double_t esdpid[15];
2396  for(Int_t i=0;i<5;i++) {
2397  esdpid[i] = esdpid0[i];
2398  esdpid[5+i] = esdpid1[i];
2399  esdpid[10+i] = esdpid2[i];
2400  }
2401  rd->SetPID(3,esdpid);
2402  return rd;
2403 }
2404 //----------------------------------------------------------------------------
2406  TObjArray *fourTrackArray,AliVEvent *event,
2407  AliAODVertex *secVert,
2408  const AliAODVertex *vertexp1n1,
2409  const AliAODVertex *vertexp1n1p2,
2410  Double_t dcap1n1,Double_t dcap1n2,
2411  Double_t dcap2n1,Double_t dcap2n2,
2412  Bool_t &ok4Prong)
2413 {
2416  // G.E.Bruno, R.Romita
2417  // AliCodeTimerAuto("",0);
2418 
2419  UInt_t ntref=TProcessID::GetObjectCount();
2420  if(ntref>16776216){// This number is 2^24-1000. The maximum number of TRef for a given TProcesssID is 2^24=16777216.
2421  AliFatal(Form("Number of TRef created (=%d), close to limit (16777216)",ntref));
2422  }
2423 
2424  ok4Prong=kFALSE;
2425  if(!secVert || !vertexp1n1 || !vertexp1n1p2) return 0x0;
2426 
2427  Double_t px[4],py[4],pz[4],d0[4],d0err[4];//d0z[3];
2428 
2429  AliESDtrack *postrack1 = (AliESDtrack*)fourTrackArray->UncheckedAt(0);
2430  AliESDtrack *negtrack1 = (AliESDtrack*)fourTrackArray->UncheckedAt(1);
2431  AliESDtrack *postrack2 = (AliESDtrack*)fourTrackArray->UncheckedAt(2);
2432  AliESDtrack *negtrack2 = (AliESDtrack*)fourTrackArray->UncheckedAt(3);
2433 
2434  postrack1->PropagateToDCA(secVert,fBzkG,kVeryBig);
2435  negtrack1->PropagateToDCA(secVert,fBzkG,kVeryBig);
2436  postrack2->PropagateToDCA(secVert,fBzkG,kVeryBig);
2437  negtrack2->PropagateToDCA(secVert,fBzkG,kVeryBig);
2438 
2439  Double_t momentum[3];
2440  postrack1->GetPxPyPz(momentum);
2441  px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2];
2442  negtrack1->GetPxPyPz(momentum);
2443  px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2];
2444  postrack2->GetPxPyPz(momentum);
2445  px[2] = momentum[0]; py[2] = momentum[1]; pz[2] = momentum[2];
2446  negtrack2->GetPxPyPz(momentum);
2447  px[3] = momentum[0]; py[3] = momentum[1]; pz[3] = momentum[2];
2448 
2449  // invariant mass cut for rho or D0 (try to improve coding here..)
2450  Bool_t okMassCut=kFALSE;
2451  if(fMassCutBeforeVertexing) okMassCut=kTRUE; // mass cut already done and passed
2452  if(!okMassCut && !(fCutsD0toKpipipi->GetUsePID())) { //no PID, to be implemented with PID
2453  if(SelectInvMassAndPt4prong(px,py,pz)) okMassCut=kTRUE;
2454  }
2455  if(!okMassCut) {
2456  //if(fDebug) printf(" candidate didn't pass mass cut\n");
2457  //printf(" candidate didn't pass mass cut\n");
2458  return 0x0;
2459  }
2460 
2461  // primary vertex to be used by this candidate
2462  AliAODVertex *primVertexAOD = PrimaryVertex(fourTrackArray,event);
2463  if(!primVertexAOD) return 0x0;
2464 
2465  Double_t d0z0[2],covd0z0[3];
2466  postrack1->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2467  d0[0]=d0z0[0];
2468  d0err[0] = TMath::Sqrt(covd0z0[0]);
2469  negtrack1->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2470  d0[1]=d0z0[0];
2471  d0err[1] = TMath::Sqrt(covd0z0[0]);
2472  postrack2->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2473  d0[2]=d0z0[0];
2474  d0err[2] = TMath::Sqrt(covd0z0[0]);
2475  negtrack2->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
2476  d0[3]=d0z0[0];
2477  d0err[3] = TMath::Sqrt(covd0z0[0]);
2478 
2479 
2480  // create the object AliAODRecoDecayHF4Prong
2481  Double_t pos[3]; primVertexAOD->GetXYZ(pos);
2482  Double_t dca[6]={dcap1n1,0.,dcap1n2,dcap2n1,0.,dcap2n2};
2483  Double_t dist12=TMath::Sqrt((vertexp1n1->GetX()-pos[0])*(vertexp1n1->GetX()-pos[0])+(vertexp1n1->GetY()-pos[1])*(vertexp1n1->GetY()-pos[1])+(vertexp1n1->GetZ()-pos[2])*(vertexp1n1->GetZ()-pos[2]));
2484  Double_t dist3=TMath::Sqrt((vertexp1n1p2->GetX()-pos[0])*(vertexp1n1p2->GetX()-pos[0])+(vertexp1n1p2->GetY()-pos[1])*(vertexp1n1p2->GetY()-pos[1])+(vertexp1n1p2->GetZ()-pos[2])*(vertexp1n1p2->GetZ()-pos[2]));
2485  Double_t dist4=TMath::Sqrt((secVert->GetX()-pos[0])*(secVert->GetX()-pos[0])+(secVert->GetY()-pos[1])*(secVert->GetY()-pos[1])+(secVert->GetZ()-pos[2])*(secVert->GetZ()-pos[2]));
2486  Short_t charge=0;
2487  AliAODRecoDecayHF4Prong *the4Prong = new AliAODRecoDecayHF4Prong(secVert,px,py,pz,d0,d0err,dca,dist12,dist3,dist4,charge);
2488  the4Prong->SetOwnPrimaryVtx(primVertexAOD);
2489  UShort_t id[4]={(UShort_t)postrack1->GetID(),(UShort_t)negtrack1->GetID(),(UShort_t)postrack2->GetID(),(UShort_t)negtrack2->GetID()};
2490  the4Prong->SetProngIDs(4,id);
2491 
2492  delete primVertexAOD; primVertexAOD=NULL;
2493 
2494  ok4Prong=(Bool_t)fCutsD0toKpipipi->IsSelected(the4Prong,AliRDHFCuts::kCandidate);
2495 
2496 
2498  the4Prong->UnsetOwnPrimaryVtx();
2499  }
2500 
2501 
2502  // get PID info from ESD
2503  Double_t esdpid0[5]={0.,0.,0.,0.,0.};
2504  if(postrack1->GetStatus()&AliESDtrack::kESDpid) postrack1->GetESDpid(esdpid0);
2505  Double_t esdpid1[5]={0.,0.,0.,0.,0.};
2506  if(negtrack1->GetStatus()&AliESDtrack::kESDpid) negtrack1->GetESDpid(esdpid1);
2507  Double_t esdpid2[5]={0.,0.,0.,0.,0.};
2508  if(postrack2->GetStatus()&AliESDtrack::kESDpid) postrack2->GetESDpid(esdpid2);
2509  Double_t esdpid3[5]={0.,0.,0.,0.,0.};
2510  if(negtrack2->GetStatus()&AliESDtrack::kESDpid) negtrack2->GetESDpid(esdpid3);
2511 
2512  Double_t esdpid[20];
2513  for(Int_t i=0;i<5;i++) {
2514  esdpid[i] = esdpid0[i];
2515  esdpid[5+i] = esdpid1[i];
2516  esdpid[10+i] = esdpid2[i];
2517  esdpid[15+i] = esdpid3[i];
2518  }
2519  the4Prong->SetPID(4,esdpid);
2520 
2521  return the4Prong;
2522 }
2523 //----------------------------------------------------------------------------------
2525  //assign and save in fAODMap the index of the AliAODTrack track
2526  //ordering them on the basis of selected criteria
2527 
2528  fAODMapSize = 100000;
2529  fAODMap = new Int_t[fAODMapSize];
2530  AliAODTrack *track=0;
2531  memset(fAODMap,0,sizeof(Int_t)*fAODMapSize);
2532  for(Int_t i=0; i<aod->GetNumberOfTracks(); i++) {
2533  track = dynamic_cast<AliAODTrack*>(aod->GetTrack(i));
2534  if(!track) AliFatal("Not a standard AOD");
2535  // skip pure ITS SA tracks
2536  if(track->GetStatus()&AliESDtrack::kITSpureSA) continue;
2537 
2538  // skip tracks without ITS
2539  if(!(track->GetStatus()&AliESDtrack::kITSin)) continue;
2540 
2541  // TEMPORARY: check that the cov matrix is there
2542  Double_t covtest[21];
2543  if(!track->GetCovarianceXYZPxPyPz(covtest)) continue;
2544  //
2545 
2546  Int_t ind = (Int_t)track->GetID();
2547  if (ind>-1 && ind < fAODMapSize) fAODMap[ind] = i;
2548  }
2549  return;
2550 }
2551 //-----------------------------------------------------------------------------
2552 AliAODVertex* AliAnalysisVertexingHF::PrimaryVertex(const TObjArray *trkArray,
2553  AliVEvent *event) const
2554 {
2556  //AliCodeTimerAuto("",0);
2557 
2558  AliESDVertex *vertexESD = 0;
2559  AliAODVertex *vertexAOD = 0;
2560 
2561 
2563  // primary vertex from the input event
2564 
2565  vertexESD = new AliESDVertex(*fV1);
2566 
2567  } else {
2568  // primary vertex specific to this candidate
2569 
2570  Int_t nTrks = trkArray->GetEntriesFast();
2571  AliVertexerTracks *vertexer = new AliVertexerTracks(event->GetMagneticField());
2572 
2574  // recalculating the vertex
2575 
2576  if(strstr(fV1->GetTitle(),"VertexerTracksWithConstraint")) {
2577  Float_t diamondcovxy[3];
2578  event->GetDiamondCovXY(diamondcovxy);
2579  Double_t pos[3]={event->GetDiamondX(),event->GetDiamondY(),0.};
2580  Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
2581  AliESDVertex *diamond = new AliESDVertex(pos,cov,1.,1);
2582  vertexer->SetVtxStart(diamond);
2583  delete diamond; diamond=NULL;
2584  if(strstr(fV1->GetTitle(),"VertexerTracksWithConstraintOnlyFitter"))
2585  vertexer->SetOnlyFitter();
2586  }
2587  Int_t skipped[1000];
2588  Int_t nTrksToSkip=0,id;
2589  AliExternalTrackParam *t = 0;
2590  for(Int_t i=0; i<nTrks; i++) {
2591  t = (AliExternalTrackParam*)trkArray->UncheckedAt(i);
2592  id = (Int_t)t->GetID();
2593  if(id<0) continue;
2594  skipped[nTrksToSkip++] = id;
2595  }
2596  // TEMPORARY FIX
2597  // For AOD, skip also tracks without covariance matrix
2598  if(fInputAOD) {
2599  Double_t covtest[21];
2600  for(Int_t j=0; j<event->GetNumberOfTracks(); j++) {
2601  AliVTrack *vtrack = (AliVTrack*)event->GetTrack(j);
2602  if(!vtrack->GetCovarianceXYZPxPyPz(covtest)) {
2603  id = (Int_t)vtrack->GetID();
2604  if(id<0) continue;
2605  skipped[nTrksToSkip++] = id;
2606  }
2607  }
2608  }
2609  for(Int_t ijk=nTrksToSkip; ijk<1000; ijk++) skipped[ijk]=-1;
2610  //
2611  vertexer->SetSkipTracks(nTrksToSkip,skipped);
2612  vertexESD = (AliESDVertex*)vertexer->FindPrimaryVertex(event);
2613 
2614  } else if(fRmTrksFromPrimVtx && nTrks>0) {
2615  // removing the prongs tracks
2616 
2617  TObjArray rmArray(nTrks);
2618  UShort_t *rmId = new UShort_t[nTrks];
2619  AliESDtrack *esdTrack = 0;
2620  AliESDtrack *t = 0;
2621  for(Int_t i=0; i<nTrks; i++) {
2622  t = (AliESDtrack*)trkArray->UncheckedAt(i);
2623  esdTrack = new AliESDtrack(*t);
2624  rmArray.AddLast(esdTrack);
2625  if(esdTrack->GetID()>=0) {
2626  rmId[i]=(UShort_t)esdTrack->GetID();
2627  } else {
2628  rmId[i]=9999;
2629  }
2630  }
2631  Float_t diamondxy[2]={static_cast<Float_t>(event->GetDiamondX()),static_cast<Float_t>(event->GetDiamondY())};
2632  vertexESD = vertexer->RemoveTracksFromVertex(fV1,&rmArray,rmId,diamondxy);
2633  delete [] rmId; rmId=NULL;
2634  rmArray.Delete();
2635 
2636  }
2637 
2638  if(!vertexESD) return vertexAOD;
2639  if(vertexESD->GetNContributors()<=0) {
2640  //AliDebug(2,"vertexing failed");
2641  delete vertexESD; vertexESD=NULL;
2642  return vertexAOD;
2643  }
2644 
2645  delete vertexer; vertexer=NULL;
2646 
2647  }
2648 
2649  // convert to AliAODVertex
2650  Double_t pos[3],cov[6],chi2perNDF;
2651  vertexESD->GetXYZ(pos); // position
2652  vertexESD->GetCovMatrix(cov); //covariance matrix
2653  chi2perNDF = vertexESD->GetChi2toNDF();
2654  delete vertexESD; vertexESD=NULL;
2655 
2656  vertexAOD = new AliAODVertex(pos,cov,chi2perNDF);
2657 
2658  return vertexAOD;
2659 }
2660 //-----------------------------------------------------------------------------
2663 
2664  //printf("Preselections:\n");
2665  // fTrackFilter->Dump();
2666  if(fSecVtxWithKF) {
2667  printf("Secondary vertex with Kalman filter package (AliKFParticle)\n");
2668  } else {
2669  printf("Secondary vertex with AliVertexerTracks\n");
2670  }
2671  if(fRecoPrimVtxSkippingTrks) printf("RecoPrimVtxSkippingTrks\n");
2672  if(fRmTrksFromPrimVtx) printf("RmTrksFromPrimVtx\n");
2673  if(fD0toKpi) {
2674  printf("Reconstruct D0->Kpi candidates with cuts:\n");
2676  }
2677  if(fDstar) {
2678  printf("Reconstruct D*->D0pi candidates with cuts:\n");
2679  if(fFindVertexForDstar) {
2680  printf(" Reconstruct a secondary vertex for the D*\n");
2681  } else {
2682  printf(" Assume the D* comes from the primary vertex\n");
2683  }
2685  }
2686  if(fJPSItoEle) {
2687  printf("Reconstruct J/psi from B candidates with cuts:\n");
2689  }
2690  if(f3Prong) {
2691  printf("Reconstruct 3 prong candidates.\n");
2692  printf(" D+->Kpipi cuts:\n");
2694  printf(" Ds->KKpi cuts:\n");
2696  printf(" Lc->pKpi cuts:\n");
2698  }
2699  if(f4Prong) {
2700  printf("Reconstruct 4 prong candidates.\n");
2701  printf(" D0->Kpipipi cuts:\n");
2703  }
2704  if(fCascades) {
2705  printf("Reconstruct cascades candidates formed with v0s.\n");
2706  printf(" Lc -> k0s P & Lc -> L Pi cuts:\n");
2708  }
2709 
2710  return;
2711 }
2712 //-----------------------------------------------------------------------------
2713 AliAODVertex* AliAnalysisVertexingHF::ReconstructSecondaryVertex(TObjArray *trkArray,
2714  Double_t &dispersion,Bool_t useTRefArray) const
2715 {
2717  //AliCodeTimerAuto("",0);
2718 
2719  AliESDVertex *vertexESD = 0;
2720  AliAODVertex *vertexAOD = 0;
2721 
2722  if(!fSecVtxWithKF) { // AliVertexerTracks
2723 
2724  fVertexerTracks->SetVtxStart(fV1);
2725  vertexESD = (AliESDVertex*)fVertexerTracks->VertexForSelectedESDTracks(trkArray);
2726 
2727  if(!vertexESD) return vertexAOD;
2728 
2729  if(vertexESD->GetNContributors()!=trkArray->GetEntriesFast()) {
2730  //AliDebug(2,"vertexing failed");
2731  delete vertexESD; vertexESD=NULL;
2732  return vertexAOD;
2733  }
2734 
2735  Double_t vertRadius2=vertexESD->GetX()*vertexESD->GetX()+vertexESD->GetY()*vertexESD->GetY();
2736  if(vertRadius2>8.){
2737  // vertex outside beam pipe, reject candidate to avoid propagation through material
2738  delete vertexESD; vertexESD=NULL;
2739  return vertexAOD;
2740  }
2741 
2742  } else { // Kalman Filter vertexer (AliKFParticle)
2743 
2744  AliKFParticle::SetField(fBzkG);
2745 
2746  AliKFVertex vertexKF;
2747 
2748  Int_t nTrks = trkArray->GetEntriesFast();
2749  for(Int_t i=0; i<nTrks; i++) {
2750  AliESDtrack *esdTrack = (AliESDtrack*)trkArray->At(i);
2751  AliKFParticle daughterKF(*esdTrack,211);
2752  vertexKF.AddDaughter(daughterKF);
2753  }
2754  vertexESD = new AliESDVertex(vertexKF.Parameters(),
2755  vertexKF.CovarianceMatrix(),
2756  vertexKF.GetChi2(),
2757  vertexKF.GetNContributors());
2758 
2759  }
2760 
2761  // convert to AliAODVertex
2762  Double_t pos[3],cov[6],chi2perNDF;
2763  vertexESD->GetXYZ(pos); // position
2764  vertexESD->GetCovMatrix(cov); //covariance matrix
2765  chi2perNDF = vertexESD->GetChi2toNDF();
2766  dispersion = vertexESD->GetDispersion();
2767  delete vertexESD; vertexESD=NULL;
2768 
2769  Int_t nprongs= (useTRefArray ? 0 : trkArray->GetEntriesFast());
2770  vertexAOD = new AliAODVertex(pos,cov,chi2perNDF,0x0,-1,AliAODVertex::kUndef,nprongs);
2771 
2772  return vertexAOD;
2773 }
2774 //-----------------------------------------------------------------------------
2777  //AliCodeTimerAuto("",0);
2778 
2779  Int_t retval=kFALSE;
2780  Double_t momentum[3];
2781  Double_t px[3],py[3],pz[3];
2782  for(Int_t iTrack=0; iTrack<3; iTrack++){
2783  AliESDtrack *track = (AliESDtrack*)trkArray->UncheckedAt(iTrack);
2784  track->GetPxPyPz(momentum);
2785  px[iTrack] = momentum[0]; py[iTrack] = momentum[1]; pz[iTrack] = momentum[2];
2786  }
2787  retval = SelectInvMassAndPt3prong(px,py,pz);
2788 
2789  return retval;
2790 }
2791 
2792 //-----------------------------------------------------------------------------
2795  //AliCodeTimerAuto("",0);
2796 
2797  Int_t retval=kFALSE;
2798  Double_t momentum[3];
2799  Double_t px[4],py[4],pz[4];
2800 
2801  for(Int_t iTrack=0; iTrack<4; iTrack++){
2802  AliESDtrack *track = (AliESDtrack*)trkArray->UncheckedAt(iTrack);
2803  track->GetPxPyPz(momentum);
2804  px[iTrack] = momentum[0]; py[iTrack] = momentum[1]; pz[iTrack] = momentum[2];
2805  }
2806 
2807  retval = SelectInvMassAndPt4prong(px,py,pz);
2808 
2809  return retval;
2810 }
2811 //-----------------------------------------------------------------------------
2814  //AliCodeTimerAuto("",0);
2815 
2816  Int_t retval=kFALSE;
2817  Double_t momentum[3];
2818  Double_t px[2],py[2],pz[2];
2819 
2820  for(Int_t iTrack=0; iTrack<2; iTrack++){
2821  AliESDtrack *track = (AliESDtrack*)trkArray->UncheckedAt(iTrack);
2822  track->GetPxPyPz(momentum);
2823  px[iTrack] = momentum[0]; py[iTrack] = momentum[1]; pz[iTrack] = momentum[2];
2824  }
2825  retval = SelectInvMassAndPtDstarD0pi(px,py,pz);
2826 
2827  return retval;
2828 }
2829 //-----------------------------------------------------------------------------
2831  Double_t *py,
2832  Double_t *pz){
2834  //AliCodeTimerAuto("",0);
2835 
2836  UInt_t pdg2[2];
2837  Int_t nprongs=2;
2838  Double_t minv2,mrange;
2839  Double_t lolim,hilim;
2840  Double_t minPt=0;
2841  Bool_t retval=kFALSE;
2842 
2843  fMassCalc2->SetPxPyPzProngs(nprongs,px,py,pz);
2844  fOKInvMassD0=kFALSE;
2845  // pt cut
2847  if(minPt>0.1)
2848  if(fMassCalc2->Pt2() < minPt*minPt) return retval;
2849  // mass cut
2850  mrange=fCutsD0toKpi->GetMassCut();
2851  lolim=fMassDzero-mrange;
2852  hilim=fMassDzero+mrange;
2853  pdg2[0]=211; pdg2[1]=321;
2854  minv2 = fMassCalc2->InvMass2(nprongs,pdg2);
2855  if(minv2>lolim*lolim && minv2<hilim*hilim ){
2856  retval=kTRUE;
2857  fOKInvMassD0=kTRUE;
2858  }
2859  pdg2[0]=321; pdg2[1]=211;
2860  minv2 = fMassCalc2->InvMass2(nprongs,pdg2);
2861  if(minv2>lolim*lolim && minv2<hilim*hilim ){
2862  retval=kTRUE;
2863  fOKInvMassD0=kTRUE;
2864  }
2865  return retval;
2866 }
2867 
2868 //-----------------------------------------------------------------------------
2870  Double_t *py,
2871  Double_t *pz){
2873  //AliCodeTimerAuto("",0);
2874 
2875  UInt_t pdg2[2];
2876  Int_t nprongs=2;
2877  Double_t minv2,mrange;
2878  Double_t lolim,hilim;
2879  Double_t minPt=0;
2880  Bool_t retval=kFALSE;
2881 
2882  fMassCalc2->SetPxPyPzProngs(nprongs,px,py,pz);
2883  fOKInvMassJpsi=kFALSE;
2884  // pt cut
2886  if(minPt>0.1)
2887  if(fMassCalc2->Pt2() < minPt*minPt) return retval;
2888  // mass cut
2889  mrange=fCutsJpsitoee->GetMassCut();
2890  lolim=fMassJpsi-mrange;
2891  hilim=fMassJpsi+mrange;
2892 
2893  pdg2[0]=11; pdg2[1]=11;
2894  minv2 = fMassCalc2->InvMass2(nprongs,pdg2);
2895  if(minv2>lolim*lolim && minv2<hilim*hilim ){
2896  retval=kTRUE;
2897  fOKInvMassJpsi=kTRUE;
2898  }
2899 
2900  return retval;
2901 }
2902 //-----------------------------------------------------------------------------
2904  Double_t *py,
2905  Double_t *pz,
2906  Int_t pidLcStatus){
2908  //AliCodeTimerAuto("",0);
2909 
2910  UInt_t pdg3[3];
2911  Int_t nprongs=3;
2912  Double_t minv2,mrange;
2913  Double_t lolim,hilim;
2914  Double_t minPt=0;
2915  Bool_t retval=kFALSE;
2916 
2917 
2918  fMassCalc3->SetPxPyPzProngs(nprongs,px,py,pz);
2919  fOKInvMassDplus=kFALSE;
2920  fOKInvMassDs=kFALSE;
2921  fOKInvMassLc=kFALSE;
2922  // pt cut
2924  minPt=TMath::Min(minPt,fCutsLctopKpi->GetMinPtCandidate());
2925  if(minPt>0.1)
2926  if(fMassCalc3->Pt2() < minPt*minPt) return retval;
2927  // D+->Kpipi
2928  mrange=fCutsDplustoKpipi->GetMassCut();
2929  lolim=fMassDplus-mrange;
2930  hilim=fMassDplus+mrange;
2931  pdg3[0]=211; pdg3[1]=321; pdg3[2]=211;
2932  minv2 = fMassCalc3->InvMass2(nprongs,pdg3);
2933  if(minv2>lolim*lolim && minv2<hilim*hilim ){
2934  retval=kTRUE;
2935  fOKInvMassDplus=kTRUE;
2936  }
2937  // Ds+->KKpi
2938  mrange=fCutsDstoKKpi->GetMassCut();
2939  lolim=fMassDs-mrange;
2940  hilim=fMassDs+mrange;
2941  pdg3[0]=321; pdg3[1]=321; pdg3[2]=211;
2942  minv2 = fMassCalc3->InvMass2(nprongs,pdg3);
2943  if(minv2>lolim*lolim && minv2<hilim*hilim ){
2944  retval=kTRUE;
2945  fOKInvMassDs=kTRUE;
2946  }
2947  pdg3[0]=211; pdg3[1]=321; pdg3[2]=321;
2948  minv2 = fMassCalc3->InvMass2(nprongs,pdg3);
2949  if(minv2>lolim*lolim && minv2<hilim*hilim ){
2950  retval=kTRUE;
2951  fOKInvMassDs=kTRUE;
2952  }
2953  // Lc->pKpi
2954  mrange=fCutsLctopKpi->GetMassCut();
2955  lolim=fMassLambdaC-mrange;
2956  hilim=fMassLambdaC+mrange;
2957  if(pidLcStatus&1){
2958  pdg3[0]=2212; pdg3[1]=321; pdg3[2]=211;
2959  minv2 = fMassCalc3->InvMass2(nprongs,pdg3);
2960  if(minv2>lolim*lolim && minv2<hilim*hilim ){
2961  retval=kTRUE;
2962  fOKInvMassLc=kTRUE;
2963  }
2964  }
2965  if(pidLcStatus&2){
2966  pdg3[0]=211; pdg3[1]=321; pdg3[2]=2212;
2967  minv2 = fMassCalc3->InvMass2(nprongs,pdg3);
2968  if(minv2>lolim*lolim && minv2<hilim*hilim ){
2969  retval=kTRUE;
2970  fOKInvMassLc=kTRUE;
2971  }
2972  }
2973 
2974  return retval;
2975 }
2976 
2977 //-----------------------------------------------------------------------------
2979  Double_t *py,
2980  Double_t *pz){
2982  //AliCodeTimerAuto("",0);
2983 
2984  UInt_t pdg2[2];
2985  Int_t nprongs=2;
2986  Double_t minv2,mrange;
2987  Double_t lolim,hilim;
2988  Double_t minPt=0;
2989  Bool_t retval=kFALSE;
2990 
2991  fMassCalc2->SetPxPyPzProngs(nprongs,px,py,pz);
2992  fOKInvMassDstar=kFALSE;
2993  // pt cut
2995  if(minPt>0.1)
2996  if(fMassCalc2->Pt2() < minPt*minPt) return retval;
2997  // mass cut
2998  pdg2[0]=211; pdg2[1]=421; // in twoTrackArrayCasc we put the pion first
2999  mrange=fCutsDStartoKpipi->GetMassCut();
3000  lolim=fMassDstar-mrange;
3001  hilim=fMassDstar+mrange;
3002  minv2 = fMassCalc2->InvMass2(nprongs,pdg2);
3003  if(minv2>lolim*lolim && minv2<hilim*hilim ){
3004  retval=kTRUE;
3005  fOKInvMassDstar=kTRUE;
3006  }
3007 
3008  return retval;
3009 }
3010 
3011 //-----------------------------------------------------------------------------
3013  Double_t *py,
3014  Double_t *pz){
3016  //AliCodeTimerAuto("",0);
3017 
3018  UInt_t pdg4[4];
3019  Int_t nprongs=4;
3020  Double_t minv2,mrange;
3021  Double_t lolim,hilim;
3022  Double_t minPt=0;
3023  Bool_t retval=kFALSE;
3024 
3025  // D0->Kpipipi without PID
3026  fMassCalc4->SetPxPyPzProngs(nprongs,px,py,pz);
3027  fOKInvMassD0to4p=kFALSE;
3028  // pt cut
3030  if(minPt>0.1)
3031  if(fMassCalc4->Pt2() < minPt*minPt) return retval;
3032  // mass cut
3033  mrange=fCutsD0toKpipipi->GetMassCut();
3034  lolim=fMassDzero-mrange;
3035  hilim=fMassDzero+mrange;
3036 
3037  pdg4[0]=321; pdg4[1]=211; pdg4[2]=211; pdg4[3]=211;
3038  minv2 = fMassCalc4->InvMass2(nprongs,pdg4);
3039  if(minv2>lolim*lolim && minv2<hilim*hilim ){
3040  retval=kTRUE;
3041  fOKInvMassD0to4p=kTRUE;
3042  }
3043 
3044  pdg4[0]=211; pdg4[1]=321; pdg4[2]=211; pdg4[3]=211;
3045  minv2 = fMassCalc4->InvMass2(nprongs,pdg4);
3046  if(minv2>lolim*lolim && minv2<hilim*hilim ){
3047  retval=kTRUE;
3048  fOKInvMassD0to4p=kTRUE;
3049  }
3050 
3051  pdg4[0]=211; pdg4[1]=211; pdg4[2]=321; pdg4[3]=211;
3052  minv2 = fMassCalc4->InvMass2(nprongs,pdg4);
3053  if(minv2>lolim*lolim && minv2<hilim*hilim ){
3054  retval=kTRUE;
3055  fOKInvMassD0to4p=kTRUE;
3056  }
3057 
3058  pdg4[0]=211; pdg4[1]=211; pdg4[2]=211; pdg4[3]=321;
3059  minv2 = fMassCalc4->InvMass2(nprongs,pdg4);
3060  if(minv2>lolim*lolim && minv2<hilim*hilim ){
3061  retval=kTRUE;
3062  fOKInvMassD0to4p=kTRUE;
3063  }
3064 
3065  return retval;
3066 }
3067 //-----------------------------------------------------------------------------
3069  Double_t *py,
3070  Double_t *pz){
3072  //AliCodeTimerAuto("",0);
3073 
3074  UInt_t pdg2[2];
3075  Int_t nprongs=2;
3076  Double_t minv2,mrange;
3077  Double_t lolim,hilim;
3078  // Double_t minPt=0;
3079  Bool_t retval=kFALSE;
3080 
3081  fMassCalc2->SetPxPyPzProngs(nprongs,px,py,pz);
3082  // minPt=fCutsLctoV0->GetMinPtCandidate();
3083  fOKInvMassLctoV0=kFALSE;
3084  mrange=fCutsLctoV0->GetMassCut();
3085  lolim=fMassLambdaC-mrange;
3086  hilim=fMassLambdaC+mrange;
3087  pdg2[0]=2212;pdg2[1]=310;
3088  minv2=fMassCalc2->InvMass2(2,pdg2);
3089  if(minv2>lolim*lolim && minv2<hilim*hilim ){
3090  retval=kTRUE;
3091  fOKInvMassLctoV0=kTRUE;
3092  }
3093  pdg2[0]=211;pdg2[1]=3122;
3094  minv2=fMassCalc2->InvMass2(2,pdg2);
3095  if(minv2>lolim*lolim && minv2<hilim*hilim ){
3096  retval=kTRUE;
3097  fOKInvMassLctoV0=kTRUE;
3098  }
3099 
3100  return retval;
3101 }
3102 //-----------------------------------------------------------------------------
3104  Int_t trkEntries,
3105  TObjArray &seleTrksArray,
3106  TObjArray &tracksAtVertex,
3107  Int_t &nSeleTrks,
3108  UChar_t *seleFlags,Int_t *evtNumber)
3109 {
3115  //AliCodeTimerAuto("",0);
3116 
3117  const AliVVertex *vprimary = event->GetPrimaryVertex();
3118 
3119  if(fV1) { delete fV1; fV1=NULL; }
3120  if(fAODMap) { delete [] fAODMap; fAODMap=NULL; }
3121 
3122  Int_t nindices=0;
3123  UShort_t *indices = 0;
3124  Double_t pos[3],cov[6];
3125  const Int_t entries = event->GetNumberOfTracks();
3126  AliCentrality* cent;
3127 
3128  if(!fInputAOD) { // ESD
3129  fV1 = new AliESDVertex(*((AliESDVertex*)vprimary));
3130  cent=((AliESDEvent*)event)->GetCentrality();
3131  } else { // AOD
3132  vprimary->GetXYZ(pos);
3133  vprimary->GetCovarianceMatrix(cov);
3134  fV1 = new AliESDVertex(pos,cov,100.,100,vprimary->GetName());
3135  if(entries<=0) return;
3136  indices = new UShort_t[entries];
3137  memset(indices,0,sizeof(UShort_t)*entries);
3138  fAODMapSize = 100000;
3139  fAODMap = new Int_t[fAODMapSize];
3140  memset(fAODMap,0,sizeof(Int_t)*fAODMapSize);
3141  cent=((AliAODEvent*)event)->GetCentrality();
3142  }
3143  Float_t centperc=0.1;
3144  if(event->GetRunNumber()<244824){
3145  centperc=cent->GetCentralityPercentile("V0M");
3146  }else{
3147  AliMultSelection *multSelection = (AliMultSelection * ) event->FindListObject("MultSelection");
3148  if(multSelection){
3149  centperc=multSelection->GetMultiplicityPercentile("V0M");
3150  Int_t qual = multSelection->GetEvSelCode();
3151  if(qual == 199 ) centperc=0.1; // use central cuts as default
3152  }
3153  }
3154  Bool_t okDisplaced=kFALSE,okSoftPi=kFALSE,okFor3Prong=kFALSE;
3155  nSeleTrks=0;
3156 
3157  // transfer ITS tracks from event to arrays
3158  for(Int_t i=0; i<entries; i++) {
3159  AliVTrack *track;
3160  track = (AliVTrack*)event->GetTrack(i);
3161 
3162  // skip pure ITS SA tracks
3163  if(track->GetStatus()&AliESDtrack::kITSpureSA) continue;
3164 
3165  // skip tracks without ITS
3166  if(!(track->GetStatus()&AliESDtrack::kITSin)) continue;
3167 
3168  // skip tracks with negative ID
3169  // (these are duplicated TPC-only AOD tracks, for jet analysis...)
3170  if(track->GetID()<0) continue;
3171 
3172  // TEMPORARY: check that the cov matrix is there
3173  Double_t covtest[21];
3174  if(!track->GetCovarianceXYZPxPyPz(covtest)) continue;
3175  //
3176 
3177  if(fInputAOD) {
3178  AliAODTrack *aodt = (AliAODTrack*)track;
3179  if(aodt->GetUsedForPrimVtxFit()) {
3180  indices[nindices]=aodt->GetID(); nindices++;
3181  }
3182  Int_t ind = (Int_t)aodt->GetID();
3183  if (ind>-1 && ind < fAODMapSize) fAODMap[ind] = i;
3184  }
3185 
3186  AliESDtrack *esdt = 0;
3187 
3188  if(!fInputAOD) {
3189  esdt = (AliESDtrack*)track;
3190  } else {
3191  esdt = new AliESDtrack(track);
3192  }
3193 
3194  // single track selection
3195  okDisplaced=kFALSE; okSoftPi=kFALSE; okFor3Prong=kFALSE;
3196  if(fMixEvent && i<trkEntries){
3197  evtNumber[i]=((AliMixedEvent*)event)->EventIndex(i);
3198  const AliVVertex* eventVtx=((AliMixedEvent*)event)->GetEventVertex(i);
3199  Double_t vtxPos[3],primPos[3],primCov[6],trasl[3];
3200  eventVtx->GetXYZ(vtxPos);
3201  vprimary->GetXYZ(primPos);
3202  eventVtx->GetCovarianceMatrix(primCov);
3203  for(Int_t ind=0;ind<3;ind++){
3204  trasl[ind]=vtxPos[ind]-primPos[ind];
3205  }
3206 
3207  Bool_t isTransl=esdt->Translate(trasl,primCov);
3208  if(!isTransl) {
3209  delete esdt;
3210  esdt = NULL;
3211  continue;
3212  }
3213  }
3214 
3215  if(SingleTrkCuts(esdt,centperc,okDisplaced,okSoftPi,okFor3Prong) && nSeleTrks<trkEntries) {
3216  esdt->PropagateToDCA(fV1,fBzkG,kVeryBig);
3217  seleTrksArray.AddLast(esdt);
3218  tracksAtVertex.AddLast(new AliExternalTrackParam(*esdt));
3219  seleFlags[nSeleTrks]=0;
3220  if(okDisplaced) SETBIT(seleFlags[nSeleTrks],kBitDispl);
3221  if(okFor3Prong) SETBIT(seleFlags[nSeleTrks],kBit3Prong);
3222  if(okSoftPi) SETBIT(seleFlags[nSeleTrks],kBitSoftPi);
3223 
3224  // Check the PID
3225  SETBIT(seleFlags[nSeleTrks],kBitPionCompat);
3226  SETBIT(seleFlags[nSeleTrks],kBitKaonCompat);
3227  SETBIT(seleFlags[nSeleTrks],kBitProtonCompat);
3228  Bool_t useTPC=kTRUE;
3229  if(fUseTOFPID){
3230  Double_t nsigmatofPi= fPidResponse->NumberOfSigmasTOF(esdt,AliPID::kPion);
3231  if(nsigmatofPi>-990. && (nsigmatofPi<-fnSigmaTOFPionLow || nsigmatofPi>fnSigmaTOFPionHi)){
3232  CLRBIT(seleFlags[nSeleTrks],kBitPionCompat);
3233  }
3234  Double_t nsigmatofK= fPidResponse->NumberOfSigmasTOF(esdt,AliPID::kKaon);
3235  if(nsigmatofK>-990. && (nsigmatofK<-fnSigmaTOFKaonLow || nsigmatofK>fnSigmaTOFKaonHi)){
3236  CLRBIT(seleFlags[nSeleTrks],kBitKaonCompat);
3237  }
3238  Double_t nsigmatofP= fPidResponse->NumberOfSigmasTOF(esdt,AliPID::kProton);
3239  if(nsigmatofP>-990. && (nsigmatofP<-fnSigmaTOFProtonLow || nsigmatofP>fnSigmaTOFProtonHi)){
3240  CLRBIT(seleFlags[nSeleTrks],kBitProtonCompat);
3241  }
3242  if(fUseTPCPIDOnlyIfNoTOF && nsigmatofPi>-990.) useTPC=kFALSE;
3243  }
3244  if(useTPC && fUseTPCPID && esdt->P()<fMaxMomForTPCPid){
3245  Double_t nsigmatpcPi= fPidResponse->NumberOfSigmasTPC(esdt,AliPID::kPion);
3246  if(nsigmatpcPi>-990. && (nsigmatpcPi<-fnSigmaTPCPionLow || nsigmatpcPi>fnSigmaTPCPionHi)){
3247  CLRBIT(seleFlags[nSeleTrks],kBitPionCompat);
3248  }
3249  Double_t nsigmatpcK= fPidResponse->NumberOfSigmasTPC(esdt,AliPID::kKaon);
3250  if(nsigmatpcK>-990. && (nsigmatpcK<-fnSigmaTPCKaonLow || nsigmatpcK>fnSigmaTPCKaonHi)){
3251  CLRBIT(seleFlags[nSeleTrks],kBitKaonCompat);
3252  }
3253  Double_t nsigmatpcP= fPidResponse->NumberOfSigmasTPC(esdt,AliPID::kProton);
3254  if(nsigmatpcP>-990. && (nsigmatpcP<-fnSigmaTPCProtonLow || nsigmatpcP>fnSigmaTPCProtonHi)){
3255  CLRBIT(seleFlags[nSeleTrks],kBitProtonCompat);
3256  }
3257  }
3258  nSeleTrks++;
3259  } else {
3260  if(fInputAOD) delete esdt;
3261  esdt = NULL;
3262  continue;
3263  }
3264 
3265  } // end loop on tracks
3266 
3267  // primary vertex from AOD
3268  if(fInputAOD) {
3269  delete fV1; fV1=NULL;
3270  vprimary->GetXYZ(pos);
3271  vprimary->GetCovarianceMatrix(cov);
3272  Double_t chi2toNDF = vprimary->GetChi2perNDF();
3273  Int_t ncontr=nindices;
3274  if(!strcmp(vprimary->GetTitle(),"VertexerTracksWithContraint")) ncontr += 1;
3275  Double_t chi2=chi2toNDF*(2.*(Double_t)ncontr-3.);
3276  fV1 = new AliESDVertex(pos,cov,chi2,ncontr,vprimary->GetName());
3277  fV1->SetTitle(vprimary->GetTitle());
3278  fV1->SetIndices(nindices,indices);
3279  }
3280  if(indices) { delete [] indices; indices=NULL; }
3281 
3282 
3283  return;
3284 }
3285 //-----------------------------------------------------------------------------
3287  //
3289  //
3290  if(fUsePidTag && cuts->GetPidHF()) {
3291  Bool_t usepid=cuts->GetIsUsePID();
3292  cuts->SetUsePID(kTRUE);
3293  if(cuts->IsSelectedPID(rd))
3294  rd->SetSelectionBit(bit);
3295  cuts->SetUsePID(usepid);
3296  }
3297  return;
3298 }
3299 //-----------------------------------------------------------------------------
3300 Bool_t AliAnalysisVertexingHF::SingleTrkCuts(AliESDtrack *trk,
3301  Float_t centralityperc,
3302  Bool_t &okDisplaced,
3303  Bool_t &okSoftPi,
3304  Bool_t &okFor3Prong) const
3305 {
3307 
3308  // this is needed to store the impact parameters
3309  //AliCodeTimerAuto("",0);
3310 
3311  trk->RelateToVertex(fV1,fBzkG,kVeryBig);
3312 
3313  UInt_t selectInfo;
3314  //
3315  // Track selection, displaced tracks -- 2 prong
3316  selectInfo = 0;
3317  if(centralityperc>=0 && fMaxCentPercentileForTightCuts>=0
3319  // central PbPb events, tighter cuts
3320  selectInfo = fTrackFilter2prongCentral->IsSelected(trk);
3321  }else{
3322  // standard cuts
3323  if(fTrackFilter) {
3324  selectInfo = fTrackFilter->IsSelected(trk);
3325  }
3326  }
3327  if(selectInfo) okDisplaced=kTRUE;
3328 
3329  // Track selection, displaced tracks -- 3 prong
3330  selectInfo = 0;
3331  if(centralityperc>=0 && fMaxCentPercentileForTightCuts>=0
3333  // central PbPb events, tighter cuts
3334  selectInfo = fTrackFilter3prongCentral->IsSelected(trk);
3335  }else{
3336  // standard cuts
3337  if(fTrackFilter) {
3338  selectInfo = fTrackFilter->IsSelected(trk);
3339  }
3340  }
3341  if(selectInfo) okFor3Prong=kTRUE;
3342 
3343  // Track selection, soft pions
3344  selectInfo = 0;
3345  if(fDstar && fTrackFilterSoftPi) {
3346  selectInfo = fTrackFilterSoftPi->IsSelected(trk);
3347  }
3348  if(selectInfo) okSoftPi=kTRUE;
3349 
3350  if(okDisplaced || okSoftPi || okFor3Prong) return kTRUE;
3351 
3352  return kFALSE;
3353 }
3354 
3355 
3356 //-----------------------------------------------------------------------------
3357 AliAODv0* AliAnalysisVertexingHF::TransformESDv0toAODv0(AliESDv0 *esdV0, TObjArray *twoTrackArrayV0){
3363  //
3364  //AliCodeTimerAuto("",0);
3365  Double_t vertex[3]; esdV0->GetXYZ(vertex[0],vertex[1],vertex[2]);
3366  AliAODVertex *vertexV0 = new AliAODVertex(vertex,esdV0->GetChi2V0(),AliAODVertex::kV0,2);
3367 
3368  // create the v0 neutral track to compute the DCA to the primary vertex
3369  Double_t xyz[3], pxpypz[3];
3370  esdV0->XvYvZv(xyz);
3371  esdV0->PxPyPz(pxpypz);
3372  Double_t cv[21]; for(int i=0; i<21; i++) cv[i]=0;
3373  AliNeutralTrackParam *trackesdV0 = new AliNeutralTrackParam(xyz,pxpypz,cv,0);
3374  if(!trackesdV0) {
3375  delete vertexV0;
3376  return 0;
3377  }
3378  Double_t d0z0[2],covd0z0[3];
3379  AliAODVertex *primVertexAOD = PrimaryVertex();
3380  trackesdV0->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
3381  Double_t dcaV0ToPrimVertex = TMath::Sqrt(covd0z0[0]);
3382  // get the v0 daughters to compute their DCA to the v0 vertex and get their momentum
3383  Double_t dcaV0DaughterToPrimVertex[2];
3384  AliExternalTrackParam *posV0track = (AliExternalTrackParam*)twoTrackArrayV0->UncheckedAt(0);
3385  AliExternalTrackParam *negV0track = (AliExternalTrackParam*)twoTrackArrayV0->UncheckedAt(1);
3386  if( !posV0track || !negV0track) {
3387  if(trackesdV0) {delete trackesdV0; trackesdV0=NULL;}
3388  delete vertexV0;
3389  delete primVertexAOD;
3390  return 0;
3391  }
3392  posV0track->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
3393  // if ( covd0z0[0]<=0.) dcaV0DaughterToPrimVertex[0] = 0;
3394  // else
3395  dcaV0DaughterToPrimVertex[0] = TMath::Sqrt(covd0z0[0]);
3396  negV0track->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
3397  // if ( covd0z0[0]<=0.)dcaV0DaughterToPrimVertex[1] = 0;
3398  // else
3399  dcaV0DaughterToPrimVertex[1] = TMath::Sqrt(covd0z0[0]);
3400  Double_t dcaV0Daughters = esdV0->GetDcaV0Daughters();
3401  Double_t pmom[3],nmom[3];
3402  esdV0->GetNPxPyPz(nmom[0],nmom[1],nmom[2]);
3403  esdV0->GetPPxPyPz(pmom[0],pmom[1],pmom[2]);
3404 
3405  AliAODv0 *aodV0 = new AliAODv0(vertexV0,dcaV0Daughters,dcaV0ToPrimVertex,pmom,nmom,dcaV0DaughterToPrimVertex);
3406  aodV0->SetOnFlyStatus(esdV0->GetOnFlyStatus());
3407 
3408  delete trackesdV0;
3409  delete primVertexAOD;
3410 
3411  return aodV0;
3412 }
3413 //-----------------------------------------------------------------------------
3414 void AliAnalysisVertexingHF::SetParametersAtVertex(AliESDtrack* esdt, const AliExternalTrackParam* extpar) const{
3416  //AliCodeTimerAuto("",0);
3417 
3418  const Double_t *par=extpar->GetParameter();
3419  const Double_t *cov=extpar->GetCovariance();
3420  Double_t alpha=extpar->GetAlpha();
3421  Double_t x=extpar->GetX();
3422  esdt->Set(x,alpha,par,cov);
3423  return;
3424 }
3425 //-----------------------------------------------------------------------------
3428 
3429  fMassDzero=TDatabasePDG::Instance()->GetParticle(421)->Mass();
3430  fMassDplus=TDatabasePDG::Instance()->GetParticle(411)->Mass();
3431  fMassDs=TDatabasePDG::Instance()->GetParticle(431)->Mass();
3432  fMassLambdaC=TDatabasePDG::Instance()->GetParticle(4122)->Mass();
3433  fMassDstar=TDatabasePDG::Instance()->GetParticle(413)->Mass();
3434  fMassJpsi=TDatabasePDG::Instance()->GetParticle(443)->Mass();
3435 }
3436 //-----------------------------------------------------------------------------
3438  //
3440  //
3441 
3442 
3443  //___ Check if the V0 type from AliRDHFCutsLctoV0 is the same as the one set in the ConfigVertexingHF.C for AliAnalysisVertexingHF
3444 
3445 
3447  {
3448  printf("ERROR: V0 type doesn not match in AliAnalysisVertexingHF (%d) required in AliRDHFCutsLctoV0 (%d)\n",fV0TypeForCascadeVertex,fCutsLctoV0->GetV0Type());
3449  return kFALSE;
3450  }
3451  return kTRUE;
3452 }
3453 //-----------------------------------------------------------------------------
3454 
Int_t charge
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
AliAODRecoDecay * fMassCalc4
for 3 prong
Bool_t fUsePidTag
upper momentum limit to apply TPC PID
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
void SetSelectionBit(Int_t i)
selection map
void SetupPID(AliVEvent *event)
AliAnalysisFilter * fTrackFilterSoftPi
Track Filter for displaced vertices in PbPb central events (tighter cuts) for 3 prong (D+...
AliVertexerTracks * fVertexerTracks
Float_t GetMassCut(Int_t iPtBin=0) const
AliPIDResponse * fPidResponse
event mixing
AliAODVertex * ReconstructSecondaryVertex(TObjArray *trkArray, Double_t &dispersion, Bool_t useTRefArray=kTRUE) const
Double_t fnSigmaTOFPionHi
Low cut value on n. of sigmas for pi TOF PID.
Double_t fnSigmaTOFProtonHi
Low cut value on n. of sigmas for p TOF PID.
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
AliAODVertex * PrimaryVertex(const TObjArray *trkArray=0x0, AliVEvent *event=0x0) const
AliAnalysisFilter * fTrackFilter3prongCentral
Track Filter for displaced vertices in PbPb central events (tighter cuts) for 2 prong (D0->Kpi) ...
AliAnalysisVertexingHF & operator=(const AliAnalysisVertexingHF &source)
AliRDHFCutsDstoKKpi * fCutsDstoKKpi
D+->Kpipi cuts.
Int_t IsD0FromDStarSelected(Double_t pt, TObject *obj, Int_t selectionLevel, AliAODEvent *aod) const
Bool_t fSecVtxWithKF
z componenent of field in kG
TList * list
AliRDHFCutsLctopKpi * fCutsLctopKpi
Ds->KKpi cuts.
Bool_t SelectInvMassAndPt3prong(Double_t *px, Double_t *py, Double_t *pz, Int_t pidLcStatus=3)
Int_t fnTrksTotal
triplet fullfilling Lc inv mass selection
Double_t fnSigmaTOFProtonLow
High cut value on n. of sigmas for p TPC PID.
Bool_t fOKInvMassDplus
pair fullfilling Jpsi inv mass selection
Double_t fnSigmaTPCKaonLow
High cut value on n. of sigmas for pi TOF PID.
Bool_t FillRecoCand(AliVEvent *event, AliAODRecoDecayHF3Prong *rd3)
Bool_t fUseTPCPIDOnlyIfNoTOF
switch use/not use TOF PID
AliRDHFCutsD0toKpi * fCutsD0toKpi
Track Filter for D* soft pion.
Bool_t fMixEvent
Like-sign triplets.
virtual void DeleteRecoD()
Bool_t fUseTOFPID
switch use/not use TPC PID
Bool_t fOKInvMassD0to4p
combination fullfilling D* inv mass selection
Bool_t fOKInvMassLc
triplet fullfilling Ds inv mass selection
Float_t GetMassCut(Int_t iPtBin=0) const
Int_t * fAODMap
map between index and ID for AOD tracks
Double_t fnSigmaTPCPionLow
flag to control usage of PID tagging
virtual Int_t IsSelectedPID(AliAODRecoDecayHF *)
Definition: AliRDHFCuts.h:281
Bool_t fFindVertexForCascades
reconstruct a secondary vertex or assume it's from the primary vertex
Double_t fnSigmaTPCKaonHi
Low cut value on n. of sigmas for K TPC PID.
AliAODPidHF * GetPidHF() const
Definition: AliRDHFCuts.h:231
Bool_t fUseKaonPIDfor3Prong
PID response.
Bool_t FillRecoCasc(AliVEvent *event, AliAODRecoCascadeHF *rc, Bool_t isDStar)
AliRDHFCutsDplustoKpipi * fCutsDplustoKpipi
J/psi->ee cuts.
AliRDHFCutsD0toKpipipi * fCutsD0toKpipipi
Lc –> v0 + bachelor cuts.
Bool_t fLikeSign3prong
Like-sign pairs.
Bool_t SelectInvMassAndPtCascade(Double_t *px, Double_t *py, Double_t *pz)
AliAODRecoDecayHF4Prong * Make4Prong(TObjArray *fourTrackArray, AliVEvent *event, AliAODVertex *secVert, const AliAODVertex *vertexp1n1, const AliAODVertex *vertexp1n1p2, Double_t dcap1n1, Double_t dcap1n2, Double_t dcap2n1, Double_t dcap2n2, Bool_t &ok4Prong)
Class for cuts on AOD reconstructed D+->Kpipi.
Float_t GetDCACut(Int_t iPtBin=0) const
Double_t fnSigmaTPCProtonLow
High cut value on n. of sigmas for K TOF PID.
void SetParametersAtVertex(AliESDtrack *esdt, const AliExternalTrackParam *extpar) const
Float_t GetMassCut(Int_t iPtBin=0) const
Float_t GetDCACut(Int_t iPtBin=0) const
Double_t fMaxMomForTPCPid
use TPC PID only for tracks that without TOF
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
Int_t GetIsFilled() const
Float_t GetMassCut(Int_t iPtBin=0) const
Float_t GetMassCut(Int_t iPtBin=0) const
AliAnalysisFilter * fTrackFilter2prongCentral
Track Filter for displaced vertices.
Double_t fnSigmaTOFKaonLow
High cut value on n. of sigmas for K TPC PID.
Float_t fMaxCentPercentileForTightCuts
High cut value on n. of sigmas for p TOF PID.
Bool_t fOKInvMassDstar
triplet fullfilling Lc inv mass selection
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel, AliAODEvent *aod)
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
Double_t fBzkG
vertexer, to compute secondary vertices
Float_t GetDCACut(Int_t iPtBin=0) const
Bool_t fUsePIDforLc2V0
PID for Lambdac: 0=no, 1=proton, 2=p and pi.
void SetProngIDs(Int_t nIDs, UShort_t *id)
virtual void PrintAll() const
UShort_t GetProngID(Int_t ip) const
Float_t GetMassCut(Int_t iPtBin=0) const
Bool_t fUseKaonPIDforDs
PID for Lambdac 2 V0: 0=no, 1=proton,.
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
Bool_t fUseTPCPID
Kaon PID usage for Ds.
Bool_t fOKInvMassLctoV0
4tracks fullfilling D0 inv mass selection
void SetPrimaryVtxRef(TObject *vtx)
primary vertex
void SetOwnPrimaryVtx(const AliAODVertex *vtx)
Bool_t fFindVertexForDstar
pointer to list of cuts for output file
Double_t fnSigmaTOFKaonHi
Low cut value on n. of sigmas for K TOF PID.
AliAODRecoDecay * fMassCalc2
to go faster in PbPb
Float_t GetDCACut(Int_t iPtBin=0) const
void FindCandidates(AliVEvent *event, TClonesArray *aodVerticesHFTClArr, TClonesArray *aodD0toKpiTClArr, TClonesArray *aodJPSItoEleTClArr, TClonesArray *aodCharm3ProngTClArr, TClonesArray *aodCharm4ProngTClArr, TClonesArray *aodDstarTClArr, TClonesArray *aodCascadesTClArr, TClonesArray *aodLikeSign2ProngTClArr, TClonesArray *aodLikeSign3ProngTClArr)
void AddDaughterRefs(AliAODVertex *v, const AliVEvent *event, const TObjArray *trkArray) const
Int_t fUsePIDforLc
Kaon PID usage for 3 prongs.
Bool_t IsEventSelected(AliVEvent *event)
Bool_t fLikeSign
cascades, Lc –> v0+track
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
Bool_t SingleTrkCuts(AliESDtrack *trk, Float_t centralityperc, Bool_t &okDisplaced, Bool_t &okSoftPi, Bool_t &ok3prong) const
void SetUsePID(Bool_t flag=kTRUE)
Definition: AliRDHFCuts.h:204
virtual void PrintAll() const
Float_t GetDCACut(Int_t iPtBin=0) const
AliRDHFCutsJpsitoee * fCutsJpsitoee
D0->Kpi cuts.
Bool_t SelectInvMassAndPtJpsiee(Double_t *px, Double_t *py, Double_t *pz)
Bool_t fOKInvMassD0
for 4 prong
void Setd0errProngs(Int_t nprongs, Double_t *d0)
Bool_t fOKInvMassJpsi
pair fullfilling D0 inv mass selection
Bool_t SelectInvMassAndPtDstarD0pi(Double_t *px, Double_t *py, Double_t *pz)
Int_t fV0TypeForCascadeVertex
reconstruct a secondary vertex or assume it's from the primary vertex
virtual Int_t IsSelected(TObject *obj, Int_t selectionLevel)
Bool_t GetIsUsePID() const
Definition: AliRDHFCuts.h:254
Bool_t fMassCutBeforeVertexing
Select which V0 type we want to use for the cascas.
Int_t fAODMapSize
input from AOD (kTRUE) or ESD (kFALSE)
void FixReferences(AliAODEvent *aod)
void SetSelectionBitForPID(AliRDHFCuts *cuts, AliAODRecoDecayHF *rd, Int_t bit)
AliAODRecoCascadeHF * MakeCascade(TObjArray *twoTrackArray, AliVEvent *event, AliAODVertex *secVert, AliAODRecoDecayHF2Prong *rd2Prong, Double_t dca, Bool_t &okDstar)
void SelectTracksAndCopyVertex(const AliVEvent *event, Int_t trkEntries, TObjArray &seleTrksArray, TObjArray &tracksAtVertex, Int_t &nSeleTrks, UChar_t *seleFlags, Int_t *evtNumber)
void SetSigmaVert(Double_t sigmaVert)
Bool_t fRecoPrimVtxSkippingTrks
if kTRUE use KF vertexer, else AliVertexerTracks
AliAODRecoDecay * fMassCalc3
for 2 prong
TList * fListOfCuts
Dstar->D0pi cuts.
Float_t GetMassCut(Int_t iPtBin=0) const
Double_t fnSigmaTPCProtonHi
Low cut value on n. of sigmas for p TPC PID.
Double_t GetMinPtCandidate() const
Definition: AliRDHFCuts.h:267
Float_t GetMassCut(Int_t iPtBin=0) const
AliRDHFCutsLctoV0 * fCutsLctoV0
Lc->pKpi cuts.
Float_t GetDCACut(Int_t iPtBin=0) const
void MapAODtracks(AliVEvent *aod)
AliAODRecoDecayHF3Prong * Make3Prong(TObjArray *threeTrackArray, AliVEvent *event, AliAODVertex *secVert, Double_t dispersion, const AliAODVertex *vertexp1n1, const AliAODVertex *vertexp2n1, Double_t dcap1n1, Double_t dcap2n1, Double_t dcap1p2, Bool_t useForLc, Bool_t useForDs, Bool_t &ok3Prong)
AliRDHFCutsDStartoKpipi * fCutsDStartoKpipi
D0->Kpipipi cuts.
AliAnalysisFilter * fTrackFilter
max. centrality percentile for using tight cuts
Float_t GetDCACut(Int_t iPtBin=0) const
Bool_t SelectInvMassAndPtD0Kpi(Double_t *px, Double_t *py, Double_t *pz)
Bool_t GetUsePID(Int_t iPtBin=0) const
AliAODv0 * TransformESDv0toAODv0(AliESDv0 *esdv0, TObjArray *twoTrackArrayV0)
Bool_t SelectInvMassAndPt4prong(Double_t *px, Double_t *py, Double_t *pz)
Double_t fnSigmaTPCPionHi
Low cut value on n. of sigmas for pi TPC PID.
Bool_t fD0toKpi
primary vertex
Bool_t fOKInvMassDs
triplet fullfilling D+ inv mass selection
void SetIsFilled(Int_t filled)
void AddRefs(AliAODVertex *v, AliAODRecoDecayHF *rd, const AliVEvent *event, const TObjArray *trkArray) const
AliAODRecoDecayHF2Prong * Make2Prong(TObjArray *twoTrackArray1, AliVEvent *event, AliAODVertex *secVert, Double_t dcap1n1, Bool_t &okD0, Bool_t &okJPSI, Bool_t &okD0fromDstar, Bool_t refill=kFALSE, AliAODRecoDecayHF2Prong *rd=0x0)
Double_t fnSigmaTOFPionLow
High cut value on n. of sigmas for pi TPC PID.