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