New class for PID of HF candidates (R. Romita)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAnalysisVertexingHF.cxx
CommitLineData
6a213b59 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//----------------------------------------------------------------------------
17// Implementation of the heavy-flavour vertexing analysis class
dcb444c9 18// Candidates are stored in the AOD as objects deriving from AliAODRecoDecay.
19// To be used as a task of AliAnalysisManager by means of the interface
20// class AliAnalysisTaskSEVertexingHF.
21// An example of usage in the macro AliAnalysisTaskSEVertexingHFTest.C.
6a213b59 22//
a25935a9 23// Contact: andrea.dainese@pd.infn.it
dc963de9 24// Contributors: E.Bruna, G.E.Bruno, A.Dainese, C.Di Gliglio,
25// F.Prino, R.Romita, X.M.Zhang
6a213b59 26//----------------------------------------------------------------------------
a9b75906 27#include <Riostream.h>
b82f6d67 28#include <TFile.h>
6a213b59 29#include <TDatabasePDG.h>
30#include <TString.h>
a9b75906 31#include <TList.h>
dcb444c9 32#include "AliLog.h"
33#include "AliVEvent.h"
34#include "AliVVertex.h"
35#include "AliVTrack.h"
6a213b59 36#include "AliVertexerTracks.h"
a353d2e4 37#include "AliKFVertex.h"
dcb444c9 38#include "AliESDEvent.h"
6a213b59 39#include "AliESDVertex.h"
e18fbfa7 40#include "AliExternalTrackParam.h"
dcfa35b3 41#include "AliNeutralTrackParam.h"
dcb444c9 42#include "AliESDtrack.h"
f8fa4595 43#include "AliESDtrackCuts.h"
b056c5e3 44#include "AliAODEvent.h"
6a213b59 45#include "AliAODRecoDecay.h"
46#include "AliAODRecoDecayHF.h"
47#include "AliAODRecoDecayHF2Prong.h"
48#include "AliAODRecoDecayHF3Prong.h"
49#include "AliAODRecoDecayHF4Prong.h"
2ff20727 50#include "AliAODRecoCascadeHF.h"
f3014dc3 51#include "AliRDHFCutsD0toKpi.h"
e3d40058 52#include "AliRDHFCutsJpsitoee.h"
53#include "AliRDHFCutsDplustoKpipi.h"
54#include "AliRDHFCutsDstoKKpi.h"
55#include "AliRDHFCutsLctopKpi.h"
56#include "AliRDHFCutsD0toKpipipi.h"
f8fa4595 57#include "AliAnalysisFilter.h"
6a213b59 58#include "AliAnalysisVertexingHF.h"
a25935a9 59#include "AliMixedEvent.h"
6a213b59 60
61ClassImp(AliAnalysisVertexingHF)
62
63//----------------------------------------------------------------------------
64AliAnalysisVertexingHF::AliAnalysisVertexingHF():
dcb444c9 65fInputAOD(kFALSE),
c8ab4e4f 66fAODMapSize(0),
67fAODMap(0),
b82f6d67 68fBzkG(0.),
69fSecVtxWithKF(kFALSE),
6a213b59 70fRecoPrimVtxSkippingTrks(kFALSE),
71fRmTrksFromPrimVtx(kFALSE),
72fV1(0x0),
b82f6d67 73fD0toKpi(kTRUE),
74fJPSItoEle(kTRUE),
75f3Prong(kTRUE),
76f4Prong(kTRUE),
2ff20727 77fDstar(kTRUE),
dc963de9 78fLikeSign(kFALSE),
a25935a9 79fMixEvent(kFALSE),
2ff20727 80fTrackFilter(0x0),
eee95e18 81fTrackFilterSoftPi(0x0),
f3014dc3 82fCutsD0toKpi(0x0),
e3d40058 83fCutsJpsitoee(0x0),
84fCutsDplustoKpipi(0x0),
85fCutsDstoKKpi(0x0),
86fCutsLctopKpi(0x0),
87fCutsD0toKpipipi(0x0),
88fCutsD0fromDstar(0x0),
a9b75906 89fListOfCuts(0x0),
eee95e18 90fFindVertexForDstar(kTRUE)
6a213b59 91{
92 // Default constructor
93
94 SetD0toKpiCuts();
95 SetBtoJPSICuts();
96 SetDplusCuts();
b82f6d67 97 SetDsCuts();
98 SetLcCuts();
2ff20727 99 SetDstarCuts();
721c0b8f 100 SetD0to4ProngsCuts();
6a213b59 101}
13977a79 102//--------------------------------------------------------------------------
103AliAnalysisVertexingHF::AliAnalysisVertexingHF(const AliAnalysisVertexingHF &source) :
104TNamed(source),
dcb444c9 105fInputAOD(source.fInputAOD),
c8ab4e4f 106fAODMapSize(source.fAODMapSize),
107fAODMap(source.fAODMap),
b82f6d67 108fBzkG(source.fBzkG),
109fSecVtxWithKF(source.fSecVtxWithKF),
13977a79 110fRecoPrimVtxSkippingTrks(source.fRecoPrimVtxSkippingTrks),
111fRmTrksFromPrimVtx(source.fRmTrksFromPrimVtx),
112fV1(source.fV1),
13977a79 113fD0toKpi(source.fD0toKpi),
114fJPSItoEle(source.fJPSItoEle),
115f3Prong(source.f3Prong),
116f4Prong(source.f4Prong),
2ff20727 117fDstar(source.fDstar),
dc963de9 118fLikeSign(source.fLikeSign),
a25935a9 119fMixEvent(source.fMixEvent),
2ff20727 120fTrackFilter(source.fTrackFilter),
eee95e18 121fTrackFilterSoftPi(source.fTrackFilterSoftPi),
f3014dc3 122fCutsD0toKpi(source.fCutsD0toKpi),
e3d40058 123fCutsJpsitoee(source.fCutsJpsitoee),
124fCutsDplustoKpipi(source.fCutsDplustoKpipi),
125fCutsDstoKKpi(source.fCutsDstoKKpi),
126fCutsLctopKpi(source.fCutsLctopKpi),
127fCutsD0toKpipipi(source.fCutsD0toKpipipi),
128fCutsD0fromDstar(source.fCutsD0fromDstar),
a9b75906 129fListOfCuts(source.fListOfCuts),
eee95e18 130fFindVertexForDstar(source.fFindVertexForDstar)
13977a79 131{
132 //
133 // Copy constructor
134 //
135 for(Int_t i=0; i<9; i++) fD0toKpiCuts[i]=source.fD0toKpiCuts[i];
136 for(Int_t i=0; i<9; i++) fBtoJPSICuts[i]=source.fBtoJPSICuts[i];
137 for(Int_t i=0; i<12; i++) fDplusCuts[i]=source.fDplusCuts[i];
f7c170b4 138 for(Int_t i=0; i<14; i++) fDsCuts[i]=source.fDsCuts[i];
2ff20727 139 for(Int_t i=0; i<12; i++) fLcCuts[i]=source.fLcCuts[i];
140 for(Int_t i=0; i<5; i++) fDstarCuts[i]=source.fDstarCuts[i];
721c0b8f 141 for(Int_t i=0; i<9; i++) fD0to4ProngsCuts[i]=source.fD0to4ProngsCuts[i];
13977a79 142}
143//--------------------------------------------------------------------------
144AliAnalysisVertexingHF &AliAnalysisVertexingHF::operator=(const AliAnalysisVertexingHF &source)
145{
146 //
147 // assignment operator
148 //
149 if(&source == this) return *this;
dcb444c9 150 fInputAOD = source.fInputAOD;
c8ab4e4f 151 fAODMapSize = source.fAODMapSize;
b82f6d67 152 fBzkG = source.fBzkG;
153 fSecVtxWithKF = source.fSecVtxWithKF;
13977a79 154 fRecoPrimVtxSkippingTrks = source.fRecoPrimVtxSkippingTrks;
155 fRmTrksFromPrimVtx = source.fRmTrksFromPrimVtx;
156 fV1 = source.fV1;
13977a79 157 fD0toKpi = source.fD0toKpi;
158 fJPSItoEle = source.fJPSItoEle;
159 f3Prong = source.f3Prong;
160 f4Prong = source.f4Prong;
2ff20727 161 fDstar = source.fDstar;
dc963de9 162 fLikeSign = source.fLikeSign;
a25935a9 163 fMixEvent = source.fMixEvent;
f8fa4595 164 fTrackFilter = source.fTrackFilter;
2ff20727 165 fTrackFilterSoftPi = source.fTrackFilterSoftPi;
f3014dc3 166 fCutsD0toKpi = source.fCutsD0toKpi;
e3d40058 167 fCutsJpsitoee = source.fCutsJpsitoee;
168 fCutsDplustoKpipi = source.fCutsDplustoKpipi;
169 fCutsDstoKKpi = source.fCutsDstoKKpi;
170 fCutsLctopKpi = source.fCutsLctopKpi;
171 fCutsD0toKpipipi = source.fCutsD0toKpipipi;
172 fCutsD0fromDstar = source.fCutsD0fromDstar;
a9b75906 173 fListOfCuts = source.fListOfCuts;
eee95e18 174 fFindVertexForDstar = source.fFindVertexForDstar;
13977a79 175
176 for(Int_t i=0; i<9; i++) fD0toKpiCuts[i]=source.fD0toKpiCuts[i];
177 for(Int_t i=0; i<9; i++) fBtoJPSICuts[i]=source.fBtoJPSICuts[i];
178 for(Int_t i=0; i<12; i++) fDplusCuts[i]=source.fDplusCuts[i];
f7c170b4 179 for(Int_t i=0; i<14; i++) fDsCuts[i]=source.fDsCuts[i];
2ff20727 180 for(Int_t i=0; i<12; i++) fLcCuts[i]=source.fLcCuts[i];
181 for(Int_t i=0; i<5; i++) fDstarCuts[i]=source.fDstarCuts[i];
721c0b8f 182 for(Int_t i=0; i<9; i++) fD0to4ProngsCuts[i]=source.fD0to4ProngsCuts[i];
13977a79 183
184 return *this;
185}
6a213b59 186//----------------------------------------------------------------------------
187AliAnalysisVertexingHF::~AliAnalysisVertexingHF() {
188 // Destructor
13977a79 189 if(fV1) { delete fV1; fV1=0; }
f8fa4595 190 if(fTrackFilter) { delete fTrackFilter; fTrackFilter=0; }
2ff20727 191 if(fTrackFilterSoftPi) { delete fTrackFilterSoftPi; fTrackFilterSoftPi=0; }
f3014dc3 192 if(fCutsD0toKpi) { delete fCutsD0toKpi; fCutsD0toKpi=0; }
e3d40058 193 if(fCutsJpsitoee) { delete fCutsJpsitoee; fCutsJpsitoee=0; }
194 if(fCutsDplustoKpipi) { delete fCutsDplustoKpipi; fCutsDplustoKpipi=0; }
195 if(fCutsDstoKKpi) { delete fCutsDstoKKpi; fCutsDstoKKpi=0; }
196 if(fCutsLctopKpi) { delete fCutsLctopKpi; fCutsLctopKpi=0; }
197 if(fCutsD0toKpipipi) { delete fCutsD0toKpipipi; fCutsD0toKpipipi=0; }
198 if(fCutsD0fromDstar) { delete fCutsD0fromDstar; fCutsD0fromDstar=0; }
c8ab4e4f 199 if(fAODMap) { delete fAODMap; fAODMap=0; }
6a213b59 200}
201//----------------------------------------------------------------------------
a9b75906 202TList *AliAnalysisVertexingHF::FillListOfCuts() {
203 // Fill list of analysis cuts
204
205 TList *list = new TList();
206 list->SetOwner();
207 list->SetName("ListOfCuts");
208
209 if(fCutsD0toKpi) {
210 AliRDHFCutsD0toKpi *cutsD0toKpi = new AliRDHFCutsD0toKpi(*fCutsD0toKpi);
211 list->Add(cutsD0toKpi);
212 }
213 if(fCutsJpsitoee) {
214 AliRDHFCutsJpsitoee *cutsJpsitoee = new AliRDHFCutsJpsitoee(*fCutsJpsitoee);
215 list->Add(cutsJpsitoee);
216 }
217 if(fCutsDplustoKpipi) {
218 AliRDHFCutsDplustoKpipi *cutsDplustoKpipi = new AliRDHFCutsDplustoKpipi(*fCutsDplustoKpipi);
219 list->Add(cutsDplustoKpipi);
220 }
221 if(fCutsDstoKKpi) {
222 AliRDHFCutsDstoKKpi *cutsDstoKKpi = new AliRDHFCutsDstoKKpi(*fCutsDstoKKpi);
223 list->Add(cutsDstoKKpi);
224 }
225 if(fCutsLctopKpi) {
226 AliRDHFCutsLctopKpi *cutsLctopKpi = new AliRDHFCutsLctopKpi(*fCutsLctopKpi);
227 list->Add(cutsLctopKpi);
228 }
229 if(fCutsD0toKpipipi) {
230 AliRDHFCutsD0toKpipipi *cutsD0toKpipipi = new AliRDHFCutsD0toKpipipi(*fCutsD0toKpipipi);
231 list->Add(cutsD0toKpipipi);
232 }
233 if(fCutsD0fromDstar) {
234 AliRDHFCutsD0toKpi *cutsD0fromDstar = new AliRDHFCutsD0toKpi(*fCutsD0fromDstar);
235 list->Add(cutsD0fromDstar);
236 }
237
238 // keep a pointer to the list
239 fListOfCuts = list;
240
241 return list;
242}
243//----------------------------------------------------------------------------
dcb444c9 244void AliAnalysisVertexingHF::FindCandidates(AliVEvent *event,
245 TClonesArray *aodVerticesHFTClArr,
246 TClonesArray *aodD0toKpiTClArr,
247 TClonesArray *aodJPSItoEleTClArr,
248 TClonesArray *aodCharm3ProngTClArr,
2ff20727 249 TClonesArray *aodCharm4ProngTClArr,
dc963de9 250 TClonesArray *aodDstarTClArr,
423fb9ae 251 TClonesArray *aodLikeSign2ProngTClArr,
252 TClonesArray *aodLikeSign3ProngTClArr)
b056c5e3 253{
254 // Find heavy-flavour vertex candidates
dcb444c9 255 // Input: ESD or AOD
b056c5e3 256 // Output: AOD (additional branches added)
257
a25935a9 258 if(!fMixEvent){
259 TString evtype = event->IsA()->GetName();
260 fInputAOD = ((evtype=="AliAODEvent") ? kTRUE : kFALSE);
261 } // if we do mixing AliVEvent is a AliMixedEvent
dcb444c9 262
e15f55cb 263 if(fInputAOD) {
264 AliDebug(2,"Creating HF candidates from AOD");
265 } else {
266 AliDebug(2,"Creating HF candidates from ESD");
267 }
268
b056c5e3 269 if(!aodVerticesHFTClArr) {
270 printf("ERROR: no aodVerticesHFTClArr");
271 return;
272 }
fa925a06 273 if((fD0toKpi || fDstar) && !aodD0toKpiTClArr) {
b056c5e3 274 printf("ERROR: no aodD0toKpiTClArr");
275 return;
276 }
277 if(fJPSItoEle && !aodJPSItoEleTClArr) {
278 printf("ERROR: no aodJPSItoEleTClArr");
279 return;
280 }
281 if(f3Prong && !aodCharm3ProngTClArr) {
282 printf("ERROR: no aodCharm3ProngTClArr");
283 return;
284 }
285 if(f4Prong && !aodCharm4ProngTClArr) {
286 printf("ERROR: no aodCharm4ProngTClArr");
287 return;
288 }
2ff20727 289 if(fDstar && !aodDstarTClArr) {
290 printf("ERROR: no aodDstarTClArr");
291 return;
292 }
423fb9ae 293 if(fLikeSign && !aodLikeSign2ProngTClArr) {
294 printf("ERROR: no aodLikeSign2ProngTClArr");
295 return;
296 }
297 if(fLikeSign && f3Prong && !aodLikeSign3ProngTClArr) {
298 printf("ERROR: no aodLikeSign2ProngTClArr");
dc963de9 299 return;
300 }
b056c5e3 301
302 // delete candidates from previous event and create references
423fb9ae 303 Int_t iVerticesHF=0,iD0toKpi=0,iJPSItoEle=0,i3Prong=0,i4Prong=0,iDstar=0,iLikeSign2Prong=0,iLikeSign3Prong=0;
b056c5e3 304 aodVerticesHFTClArr->Delete();
305 iVerticesHF = aodVerticesHFTClArr->GetEntriesFast();
306 TClonesArray &verticesHFRef = *aodVerticesHFTClArr;
fa925a06 307 if(fD0toKpi || fDstar) {
b056c5e3 308 aodD0toKpiTClArr->Delete();
309 iD0toKpi = aodD0toKpiTClArr->GetEntriesFast();
310 }
311 if(fJPSItoEle) {
312 aodJPSItoEleTClArr->Delete();
313 iJPSItoEle = aodJPSItoEleTClArr->GetEntriesFast();
314 }
315 if(f3Prong) {
316 aodCharm3ProngTClArr->Delete();
317 i3Prong = aodCharm3ProngTClArr->GetEntriesFast();
318 }
319 if(f4Prong) {
320 aodCharm4ProngTClArr->Delete();
321 i4Prong = aodCharm4ProngTClArr->GetEntriesFast();
322 }
2ff20727 323 if(fDstar) {
324 aodDstarTClArr->Delete();
325 iDstar = aodDstarTClArr->GetEntriesFast();
326 }
423fb9ae 327 if(fLikeSign) {
328 aodLikeSign2ProngTClArr->Delete();
329 iLikeSign2Prong = aodLikeSign2ProngTClArr->GetEntriesFast();
330 }
331 if(fLikeSign && f3Prong) {
332 aodLikeSign3ProngTClArr->Delete();
333 iLikeSign3Prong = aodLikeSign3ProngTClArr->GetEntriesFast();
dc963de9 334 }
423fb9ae 335
336 TClonesArray &aodD0toKpiRef = *aodD0toKpiTClArr;
337 TClonesArray &aodJPSItoEleRef = *aodJPSItoEleTClArr;
338 TClonesArray &aodCharm3ProngRef = *aodCharm3ProngTClArr;
339 TClonesArray &aodCharm4ProngRef = *aodCharm4ProngTClArr;
340 TClonesArray &aodDstarRef = *aodDstarTClArr;
341 TClonesArray &aodLikeSign2ProngRef = *aodLikeSign2ProngTClArr;
342 TClonesArray &aodLikeSign3ProngRef = *aodLikeSign3ProngTClArr;
dc963de9 343
b056c5e3 344
2ff20727 345 AliAODRecoDecayHF2Prong *io2Prong = 0;
346 AliAODRecoDecayHF3Prong *io3Prong = 0;
347 AliAODRecoDecayHF4Prong *io4Prong = 0;
348 AliAODRecoCascadeHF *ioCascade = 0;
b056c5e3 349
2ff20727 350 Int_t iTrkP1,iTrkP2,iTrkN1,iTrkN2,iTrkSoftPi,trkEntries;
721c0b8f 351 Double_t xdummy,ydummy,dcap1n1,dcap1n2,dcap2n1,dcap1p2,dcan1n2,dcap2n2,dcaCasc;
b056c5e3 352 Bool_t okD0=kFALSE,okJPSI=kFALSE,ok3Prong=kFALSE,ok4Prong=kFALSE;
2ff20727 353 Bool_t okDstar=kFALSE,okD0fromDstar=kFALSE;
b056c5e3 354 AliESDtrack *postrack1 = 0;
355 AliESDtrack *postrack2 = 0;
356 AliESDtrack *negtrack1 = 0;
357 AliESDtrack *negtrack2 = 0;
2ff20727 358 AliESDtrack *trackPi = 0;
a9b75906 359 /*
b056c5e3 360 Double_t dcaMax = fD0toKpiCuts[1];
361 if(dcaMax < fBtoJPSICuts[1]) dcaMax=fBtoJPSICuts[1];
362 if(dcaMax < fDplusCuts[11]) dcaMax=fDplusCuts[11];
721c0b8f 363 if(dcaMax < fD0to4ProngsCuts[1]) dcaMax=fD0to4ProngsCuts[1];
a9b75906 364 */
365 Float_t dcaMax = fCutsD0toKpi->GetDCACut();
e3d40058 366 if(fCutsJpsitoee) dcaMax=TMath::Max(dcaMax,fCutsJpsitoee->GetDCACut());
367 if(fCutsDplustoKpipi) dcaMax=TMath::Max(dcaMax,fCutsDplustoKpipi->GetDCACut());
368 if(fCutsDstoKKpi) dcaMax=TMath::Max(dcaMax,fCutsDstoKKpi->GetDCACut());
369 if(fCutsLctopKpi) dcaMax=TMath::Max(dcaMax,fCutsLctopKpi->GetDCACut());
370 if(fCutsD0toKpipipi) dcaMax=TMath::Max(dcaMax,fCutsD0toKpipipi->GetDCACut());
a9b75906 371
dcb444c9 372 AliDebug(2,Form(" dca cut set to %f cm",dcaMax));
b056c5e3 373
b056c5e3 374
dcb444c9 375 // get Bz
376 fBzkG = (Double_t)event->GetMagneticField();
b056c5e3 377
dcb444c9 378 trkEntries = (Int_t)event->GetNumberOfTracks();
379 AliDebug(1,Form(" Number of tracks: %d",trkEntries));
357cd934 380
381 if(trkEntries<2) {
382 AliDebug(1,Form(" Not enough tracks: %d",trkEntries));
383 return;
384 }
b056c5e3 385
3cbc09cd 386 // event selection
a9b75906 387 if(!fCutsD0toKpi->IsEventSelected(event)) return;
b056c5e3 388
dcb444c9 389 // call function that applies sigle-track selection,
2ff20727 390 // for displaced tracks and soft pions (both charges) for D*,
dcb444c9 391 // and retrieves primary vertex
2ff20727 392 TObjArray seleTrksArray(trkEntries);
393 UChar_t *seleFlags = new UChar_t[trkEntries]; // bit 0: displaced, bit 1: softpi
394 Int_t nSeleTrks=0;
a25935a9 395 Int_t *evtNumber = new Int_t[trkEntries];
396 SelectTracksAndCopyVertex(event,seleTrksArray,nSeleTrks,seleFlags,evtNumber);
b056c5e3 397
2ff20727 398 AliDebug(1,Form(" Selected tracks: %d",nSeleTrks));
b056c5e3 399
2ff20727 400 TObjArray *twoTrackArray1 = new TObjArray(2);
401 TObjArray *twoTrackArray2 = new TObjArray(2);
402 TObjArray *twoTrackArrayCasc = new TObjArray(2);
403 TObjArray *threeTrackArray = new TObjArray(3);
404 TObjArray *fourTrackArray = new TObjArray(4);
b056c5e3 405
dcb444c9 406 Double_t dispersion;
423fb9ae 407 Bool_t isLikeSign2Prong=kFALSE,isLikeSign3Prong=kFALSE;
dcb444c9 408
2ff20727 409 AliAODRecoDecayHF *rd = 0;
410 AliAODRecoCascadeHF *rc = 0;
411
b056c5e3 412 // LOOP ON POSITIVE TRACKS
2ff20727 413 for(iTrkP1=0; iTrkP1<nSeleTrks; iTrkP1++) {
dc963de9 414
2ff20727 415 if(iTrkP1%1==0) AliDebug(1,Form(" 1st loop on pos: track number %d of %d",iTrkP1,nSeleTrks));
dc963de9 416
b056c5e3 417 // get track from tracks array
2ff20727 418 postrack1 = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkP1);
dc963de9 419
420 if(!TESTBIT(seleFlags[iTrkP1],kBitDispl)) continue;
421
422 if(postrack1->Charge()<0 && !fLikeSign) continue;
423
b056c5e3 424 // LOOP ON NEGATIVE TRACKS
2ff20727 425 for(iTrkN1=0; iTrkN1<nSeleTrks; iTrkN1++) {
dc963de9 426
2ff20727 427 if(iTrkN1%1==0) AliDebug(1,Form(" 1st loop on neg: track number %d of %d",iTrkN1,nSeleTrks));
dc963de9 428
429 if(iTrkN1==iTrkP1) continue;
430
b056c5e3 431 // get track from tracks array
2ff20727 432 negtrack1 = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkN1);
dc963de9 433
434 if(negtrack1->Charge()>0 && !fLikeSign) continue;
435
436 if(!TESTBIT(seleFlags[iTrkN1],kBitDispl)) continue;
437
a25935a9 438 if(fMixEvent) {
439 if(evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
440 }
441
dc963de9 442 if(postrack1->Charge()==negtrack1->Charge()) { // like-sign
423fb9ae 443 isLikeSign2Prong=kTRUE;
dc963de9 444 if(!fLikeSign) continue;
445 if(iTrkN1<iTrkP1) continue; // this is needed to avoid double-counting of like-sign
446 } else { // unlike-sign
423fb9ae 447 isLikeSign2Prong=kFALSE;
dc963de9 448 if(postrack1->Charge()<0 || negtrack1->Charge()>0) continue; // this is needed to avoid double-counting of unlike-sign
a25935a9 449 if(fMixEvent) {
450 if(evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
451 }
452
dc963de9 453 }
454
dcb444c9 455 // back to primary vertex
456 postrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
457 negtrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
dc963de9 458
b056c5e3 459 // DCA between the two tracks
460 dcap1n1 = postrack1->GetDCA(negtrack1,fBzkG,xdummy,ydummy);
461 if(dcap1n1>dcaMax) { negtrack1=0; continue; }
dc963de9 462
b056c5e3 463 // Vertexing
464 twoTrackArray1->AddAt(postrack1,0);
465 twoTrackArray1->AddAt(negtrack1,1);
dcb444c9 466 AliAODVertex *vertexp1n1 = ReconstructSecondaryVertex(twoTrackArray1,dispersion);
b056c5e3 467 if(!vertexp1n1) {
468 twoTrackArray1->Clear();
469 negtrack1=0;
470 continue;
471 }
dc963de9 472
473 // 2 prong candidate
474 if(fD0toKpi || fJPSItoEle || fDstar || fLikeSign) {
b6e66d86 475
2ff20727 476 io2Prong = Make2Prong(twoTrackArray1,event,vertexp1n1,dcap1n1,okD0,okJPSI,okD0fromDstar);
b6e66d86 477
423fb9ae 478 if((fD0toKpi && okD0) || (fJPSItoEle && okJPSI) || (isLikeSign2Prong && (okD0 || okJPSI))) {
2ff20727 479 // add the vertex and the decay to the AOD
480 AliAODVertex *v2Prong = new(verticesHFRef[iVerticesHF++])AliAODVertex(*vertexp1n1);
423fb9ae 481 if(!isLikeSign2Prong) {
dc963de9 482 if(okD0) {
483 rd = new(aodD0toKpiRef[iD0toKpi++])AliAODRecoDecayHF2Prong(*io2Prong);
484 rd->SetSecondaryVtx(v2Prong);
485 v2Prong->SetParent(rd);
a9b75906 486 AddRefs(v2Prong,rd,event,twoTrackArray1);
dc963de9 487 }
488 if(okJPSI) {
489 rd = new(aodJPSItoEleRef[iJPSItoEle++])AliAODRecoDecayHF2Prong(*io2Prong);
490 rd->SetSecondaryVtx(v2Prong);
491 if(!okD0) v2Prong->SetParent(rd); // it cannot have two mothers ...
a9b75906 492 AddRefs(v2Prong,rd,event,twoTrackArray1);
dc963de9 493 }
423fb9ae 494 } else { // isLikeSign2Prong
495 rd = new(aodLikeSign2ProngRef[iLikeSign2Prong++])AliAODRecoDecayHF2Prong(*io2Prong);
2ff20727 496 rd->SetSecondaryVtx(v2Prong);
497 v2Prong->SetParent(rd);
dcb444c9 498 }
b056c5e3 499 }
dc963de9 500 // D* candidates
423fb9ae 501 if(fDstar && okD0fromDstar && !isLikeSign2Prong) {
dc963de9 502 // write references in io2Prong
503 if(fInputAOD) {
2ff20727 504 AddDaughterRefs(vertexp1n1,event,twoTrackArray1);
dc963de9 505 } else {
506 vertexp1n1->AddDaughter(postrack1);
507 vertexp1n1->AddDaughter(negtrack1);
2ff20727 508 }
dc963de9 509 io2Prong->SetSecondaryVtx(vertexp1n1);
510 //printf("---> %d %d %d %d %d\n",vertexp1n1->GetNDaughters(),iTrkP1,iTrkN1,postrack1->Charge(),negtrack1->Charge());
2ff20727 511 // create a track from the D0
dcfa35b3 512 AliNeutralTrackParam *trackD0 = new AliNeutralTrackParam(io2Prong);
dc963de9 513
2ff20727 514 // LOOP ON TRACKS THAT PASSED THE SOFT PION CUTS
515 for(iTrkSoftPi=0; iTrkSoftPi<nSeleTrks; iTrkSoftPi++) {
dc963de9 516
2ff20727 517 if(iTrkSoftPi==iTrkP1 || iTrkSoftPi==iTrkN1) continue;
dc963de9 518
2ff20727 519 if(!TESTBIT(seleFlags[iTrkSoftPi],kBitSoftPi)) continue;
dc963de9 520
a25935a9 521 if(fMixEvent) {
522 if(evtNumber[iTrkP1]==evtNumber[iTrkSoftPi] ||
523 evtNumber[iTrkN1]==evtNumber[iTrkSoftPi] ||
524 evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
525 }
526
2ff20727 527 if(iTrkSoftPi%1==0) AliDebug(1,Form(" 1st loop on pi_s: track number %d of %d",iTrkSoftPi,nSeleTrks));
dc963de9 528
eee95e18 529 trackD0->PropagateToDCA(fV1,fBzkG,kVeryBig);
530 if(trackD0->GetSigmaY2()<0. || trackD0->GetSigmaZ2()<0.) continue; // this is insipired by the AliITStrackV2::Invariant() checks
531
2ff20727 532 // get track from tracks array
533 trackPi = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkSoftPi);
2ff20727 534 trackPi->PropagateToDCA(fV1,fBzkG,kVeryBig);
dc963de9 535
2ff20727 536 twoTrackArrayCasc->AddAt(trackPi,0);
537 twoTrackArrayCasc->AddAt(trackD0,1);
eee95e18 538
539 AliAODVertex *vertexCasc = 0;
540
541 if(fFindVertexForDstar) {
542 // DCA between the two tracks
543 dcaCasc = trackPi->GetDCA(trackD0,fBzkG,xdummy,ydummy);
544 // Vertexing
545 vertexCasc = ReconstructSecondaryVertex(twoTrackArrayCasc,dispersion,kFALSE);
546 } else {
547 // assume Dstar decays at the primary vertex
548 Double_t pos[3],cov[6],chi2perNDF;
549 fV1->GetXYZ(pos);
550 fV1->GetCovMatrix(cov);
551 chi2perNDF = fV1->GetChi2toNDF();
552 vertexCasc = new AliAODVertex(pos,cov,chi2perNDF,0x0,-1,AliAODVertex::kUndef,2);
553 dcaCasc = 0.;
554 }
2ff20727 555 if(!vertexCasc) {
556 twoTrackArrayCasc->Clear();
557 trackPi=0;
558 continue;
559 }
dc963de9 560
2ff20727 561 ioCascade = MakeCascade(twoTrackArrayCasc,event,vertexCasc,io2Prong,dcaCasc,okDstar);
562 if(okDstar) {
dc963de9 563 // add the D0 to the AOD (if not already done)
2ff20727 564 if(!okD0) {
2ff20727 565 AliAODVertex *v2Prong = new(verticesHFRef[iVerticesHF++])AliAODVertex(*vertexp1n1);
2ff20727 566 rd = new(aodD0toKpiRef[iD0toKpi++])AliAODRecoDecayHF2Prong(*io2Prong);
567 rd->SetSecondaryVtx(v2Prong);
568 v2Prong->SetParent(rd);
a9b75906 569 AddRefs(v2Prong,rd,event,twoTrackArray1);
2ff20727 570 okD0=kTRUE; // this is done to add it only once
571 }
572 // add the vertex and the cascade to the AOD
573 AliAODVertex *vCasc = new(verticesHFRef[iVerticesHF++])AliAODVertex(*vertexCasc);
2ff20727 574 rc = new(aodDstarRef[iDstar++])AliAODRecoCascadeHF(*ioCascade);
2ff20727 575 rc->SetSecondaryVtx(vCasc);
576 vCasc->SetParent(rc);
a9b75906 577 if(!fInputAOD) vCasc->AddDaughter(rd); // just to fill ref #0
578 AddRefs(vCasc,rc,event,twoTrackArrayCasc);
579 vCasc->AddDaughter(rd); // add the D0 (in ref #1)
2ff20727 580 }
581 twoTrackArrayCasc->Clear();
582 trackPi=0;
0a65d33f 583 if(ioCascade) {delete ioCascade; ioCascade=NULL;}
b6e66d86 584 delete vertexCasc; vertexCasc=NULL;
2ff20727 585 } // end loop on soft pi tracks
dc963de9 586
2ff20727 587 if(trackD0) {delete trackD0; trackD0=NULL;}
dc963de9 588
2ff20727 589 }
0a65d33f 590 if(io2Prong) {delete io2Prong; io2Prong=NULL;}
2ff20727 591 }
dc963de9 592
b056c5e3 593 twoTrackArray1->Clear();
423fb9ae 594 if( (!f3Prong && !f4Prong) ||
595 (isLikeSign2Prong && !f3Prong) ) {
b056c5e3 596 negtrack1=0;
597 delete vertexp1n1;
598 continue;
599 }
600
601
602 // 2nd LOOP ON POSITIVE TRACKS
2ff20727 603 for(iTrkP2=iTrkP1+1; iTrkP2<nSeleTrks; iTrkP2++) {
dc963de9 604
2ff20727 605 if(iTrkP2==iTrkP1 || iTrkP2==iTrkN1) continue;
dc963de9 606
2ff20727 607 if(iTrkP2%1==0) AliDebug(1,Form(" 2nd loop on pos: track number %d of %d",iTrkP2,nSeleTrks));
dc963de9 608
b056c5e3 609 // get track from tracks array
2ff20727 610 postrack2 = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkP2);
dc963de9 611
612 if(postrack2->Charge()<0) continue;
613
614 if(!TESTBIT(seleFlags[iTrkP2],kBitDispl)) continue;
615
a25935a9 616 if(fMixEvent) {
617 if(evtNumber[iTrkP1]==evtNumber[iTrkP2] ||
618 evtNumber[iTrkN1]==evtNumber[iTrkP2] ||
619 evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
620 }
621
423fb9ae 622 if(isLikeSign2Prong) { // like-sign pair -> have to build only like-sign triplet
623 if(postrack1->Charge()>0) { // ok: like-sign triplet (+++)
624 isLikeSign3Prong=kTRUE;
625 } else { // not ok
626 continue;
627 }
628 } else { // normal triplet (+-+)
629 isLikeSign3Prong=kFALSE;
a25935a9 630 if(fMixEvent) {
631 if(evtNumber[iTrkP1]==evtNumber[iTrkP2] ||
632 evtNumber[iTrkN1]==evtNumber[iTrkP2] ||
633 evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
634 }
423fb9ae 635 }
636
dcb444c9 637 // back to primary vertex
638 postrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
639 postrack2->PropagateToDCA(fV1,fBzkG,kVeryBig);
640 negtrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
2ff20727 641 //printf("********** %d %d %d\n",postrack1->GetID(),postrack2->GetID(),negtrack1->GetID());
dc963de9 642
b056c5e3 643 dcap2n1 = postrack2->GetDCA(negtrack1,fBzkG,xdummy,ydummy);
644 if(dcap2n1>dcaMax) { postrack2=0; continue; }
645 dcap1p2 = postrack2->GetDCA(postrack1,fBzkG,xdummy,ydummy);
646 if(dcap1p2>dcaMax) { postrack2=0; continue; }
647
648 // Vertexing
649 twoTrackArray2->AddAt(postrack2,0);
650 twoTrackArray2->AddAt(negtrack1,1);
dcb444c9 651 AliAODVertex *vertexp2n1 = ReconstructSecondaryVertex(twoTrackArray2,dispersion);
b056c5e3 652 if(!vertexp2n1) {
653 twoTrackArray2->Clear();
654 postrack2=0;
423fb9ae 655 continue;
b056c5e3 656 }
dc963de9 657
658 // 3 prong candidates
b056c5e3 659 if(f3Prong) {
423fb9ae 660 if(postrack2->Charge()>0) {
661 threeTrackArray->AddAt(postrack1,0);
662 threeTrackArray->AddAt(negtrack1,1);
663 threeTrackArray->AddAt(postrack2,2);
664 } else {
665 threeTrackArray->AddAt(negtrack1,0);
666 threeTrackArray->AddAt(postrack1,1);
667 threeTrackArray->AddAt(postrack2,2);
668 }
dc963de9 669
dcb444c9 670 AliAODVertex* secVert3PrAOD = ReconstructSecondaryVertex(threeTrackArray,dispersion);
671 io3Prong = Make3Prong(threeTrackArray,event,secVert3PrAOD,dispersion,vertexp1n1,vertexp2n1,dcap1n1,dcap2n1,dcap1p2,ok3Prong);
b056c5e3 672 if(ok3Prong) {
2ff20727 673 AliAODVertex *v3Prong = new(verticesHFRef[iVerticesHF++])AliAODVertex(*secVert3PrAOD);
423fb9ae 674 if(!isLikeSign3Prong) {
675 rd = new(aodCharm3ProngRef[i3Prong++])AliAODRecoDecayHF3Prong(*io3Prong);
676 rd->SetSecondaryVtx(v3Prong);
677 v3Prong->SetParent(rd);
a9b75906 678 AddRefs(v3Prong,rd,event,threeTrackArray);
423fb9ae 679 } else { // isLikeSign3Prong
680 rd = new(aodLikeSign3ProngRef[iLikeSign3Prong++])AliAODRecoDecayHF3Prong(*io3Prong);
681 rd->SetSecondaryVtx(v3Prong);
682 v3Prong->SetParent(rd);
a9b75906 683 AddRefs(v3Prong,rd,event,threeTrackArray);
423fb9ae 684 }
b056c5e3 685 }
b6e66d86 686 if(io3Prong) {delete io3Prong; io3Prong=NULL;}
687 if(secVert3PrAOD) {delete secVert3PrAOD; secVert3PrAOD=NULL;}
b056c5e3 688 }
dc963de9 689
690 // 4 prong candidates
423fb9ae 691 if(f4Prong
692 // don't make 4 prong with like-sign pairs and triplets
693 && !isLikeSign2Prong && !isLikeSign3Prong
fae0f82d 694 // track-to-track dca cuts already now
695 && dcap1n1 < fD0to4ProngsCuts[1]
696 && dcap2n1 < fD0to4ProngsCuts[1]) {
dc963de9 697
721c0b8f 698 // back to primary vertex
699 postrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
700 postrack2->PropagateToDCA(fV1,fBzkG,kVeryBig);
701 negtrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
702 // Vertexing for these 3 (can be taken from above?)
703 threeTrackArray->AddAt(postrack1,0);
704 threeTrackArray->AddAt(negtrack1,1);
705 threeTrackArray->AddAt(postrack2,2);
706 AliAODVertex* vertexp1n1p2 = ReconstructSecondaryVertex(threeTrackArray,dispersion);
707
b056c5e3 708 // 3rd LOOP ON NEGATIVE TRACKS (for 4 prong)
2ff20727 709 for(iTrkN2=iTrkN1+1; iTrkN2<nSeleTrks; iTrkN2++) {
dc963de9 710
2ff20727 711 if(iTrkN2==iTrkP1 || iTrkN2==iTrkP2 || iTrkN2==iTrkN1) continue;
dc963de9 712
2ff20727 713 if(iTrkN2%1==0) AliDebug(1,Form(" 3rd loop on neg: track number %d of %d",iTrkN2,nSeleTrks));
dc963de9 714
b056c5e3 715 // get track from tracks array
2ff20727 716 negtrack2 = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkN2);
dc963de9 717
718 if(negtrack2->Charge()>0) continue;
719
720 if(!TESTBIT(seleFlags[iTrkN2],kBitDispl)) continue;
a25935a9 721 if(fMixEvent){
722 if(evtNumber[iTrkP1]==evtNumber[iTrkN2] ||
723 evtNumber[iTrkN1]==evtNumber[iTrkN2] ||
724 evtNumber[iTrkP2]==evtNumber[iTrkN2] ||
725 evtNumber[iTrkP1]==evtNumber[iTrkN1] ||
726 evtNumber[iTrkP1]==evtNumber[iTrkP2] ||
727 evtNumber[iTrkN1]==evtNumber[iTrkP2]) continue;
728 }
dc963de9 729
dcb444c9 730 // back to primary vertex
731 postrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
732 postrack2->PropagateToDCA(fV1,fBzkG,kVeryBig);
733 negtrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
734 negtrack2->PropagateToDCA(fV1,fBzkG,kVeryBig);
b056c5e3 735 dcap1n2 = postrack1->GetDCA(negtrack2,fBzkG,xdummy,ydummy);
fae0f82d 736 if(dcap1n2 > fD0to4ProngsCuts[1]) { negtrack2=0; continue; }
721c0b8f 737 dcap2n2 = postrack2->GetDCA(negtrack2,fBzkG,xdummy,ydummy);
fae0f82d 738 if(dcap2n2 > fD0to4ProngsCuts[1]) { negtrack2=0; continue; }
dc963de9 739
b056c5e3 740 // Vertexing
741 fourTrackArray->AddAt(postrack1,0);
742 fourTrackArray->AddAt(negtrack1,1);
743 fourTrackArray->AddAt(postrack2,2);
744 fourTrackArray->AddAt(negtrack2,3);
dc963de9 745
dcb444c9 746 AliAODVertex* secVert4PrAOD = ReconstructSecondaryVertex(fourTrackArray,dispersion);
721c0b8f 747 io4Prong = Make4Prong(fourTrackArray,event,secVert4PrAOD,vertexp1n1,vertexp1n1p2,dcap1n1,dcap1n2,dcap2n1,dcap2n2,ok4Prong);
b056c5e3 748 if(ok4Prong) {
2ff20727 749 AliAODVertex *v4Prong = new(verticesHFRef[iVerticesHF++])AliAODVertex(*secVert4PrAOD);
2ff20727 750 rd = new(aodCharm4ProngRef[i4Prong++])AliAODRecoDecayHF4Prong(*io4Prong);
751 rd->SetSecondaryVtx(v4Prong);
752 v4Prong->SetParent(rd);
a9b75906 753 AddRefs(v4Prong,rd,event,fourTrackArray);
b056c5e3 754 }
dc963de9 755
b6e66d86 756 if(io4Prong) {delete io4Prong; io4Prong=NULL;}
757 if(secVert4PrAOD) {delete secVert4PrAOD; secVert4PrAOD=NULL;}
b056c5e3 758 fourTrackArray->Clear();
759 negtrack2 = 0;
dc963de9 760
b056c5e3 761 } // end loop on negative tracks
dc963de9 762
721c0b8f 763 threeTrackArray->Clear();
764 delete vertexp1n1p2;
765
b056c5e3 766 }
dc963de9 767
b056c5e3 768 postrack2 = 0;
769 delete vertexp2n1;
dc963de9 770
b056c5e3 771 } // end 2nd loop on positive tracks
dc963de9 772
b056c5e3 773 twoTrackArray2->Clear();
774
423fb9ae 775 // 2nd LOOP ON NEGATIVE TRACKS (for 3 prong -+-)
2ff20727 776 for(iTrkN2=iTrkN1+1; iTrkN2<nSeleTrks; iTrkN2++) {
dc963de9 777
2ff20727 778 if(iTrkN2==iTrkP1 || iTrkN2==iTrkP2 || iTrkN2==iTrkN1) continue;
dc963de9 779
2ff20727 780 if(iTrkN2%1==0) AliDebug(1,Form(" 2nd loop on neg: track number %d of %d",iTrkN2,nSeleTrks));
dc963de9 781
b056c5e3 782 // get track from tracks array
2ff20727 783 negtrack2 = (AliESDtrack*)seleTrksArray.UncheckedAt(iTrkN2);
dc963de9 784
785 if(negtrack2->Charge()>0) continue;
786
787 if(!TESTBIT(seleFlags[iTrkN2],kBitDispl)) continue;
788
a25935a9 789 if(fMixEvent) {
790 if(evtNumber[iTrkP1]==evtNumber[iTrkN2] ||
791 evtNumber[iTrkN1]==evtNumber[iTrkN2] ||
792 evtNumber[iTrkP1]==evtNumber[iTrkN1]) continue;
793 }
794
423fb9ae 795 if(isLikeSign2Prong) { // like-sign pair -> have to build only like-sign triplet
796 if(postrack1->Charge()<0) { // ok: like-sign triplet (---)
797 isLikeSign3Prong=kTRUE;
798 } else { // not ok
799 continue;
800 }
801 } else { // normal triplet (-+-)
a25935a9 802 isLikeSign3Prong=kFALSE;
423fb9ae 803 }
804
dcb444c9 805 // back to primary vertex
806 postrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
807 negtrack1->PropagateToDCA(fV1,fBzkG,kVeryBig);
808 negtrack2->PropagateToDCA(fV1,fBzkG,kVeryBig);
2ff20727 809 //printf("********** %d %d %d\n",postrack1->GetID(),negtrack1->GetID(),negtrack2->GetID());
dc963de9 810
b056c5e3 811 dcap1n2 = postrack1->GetDCA(negtrack2,fBzkG,xdummy,ydummy);
812 if(dcap1n2>dcaMax) { negtrack2=0; continue; }
813 dcan1n2 = negtrack1->GetDCA(negtrack2,fBzkG,xdummy,ydummy);
814 if(dcan1n2>dcaMax) { negtrack2=0; continue; }
815
816 // Vertexing
817 twoTrackArray2->AddAt(postrack1,0);
818 twoTrackArray2->AddAt(negtrack2,1);
dcb444c9 819
820 AliAODVertex *vertexp1n2 = ReconstructSecondaryVertex(twoTrackArray2,dispersion);
b056c5e3 821 if(!vertexp1n2) {
822 twoTrackArray2->Clear();
823 negtrack2=0;
824 continue;
825 }
dc963de9 826
b056c5e3 827 if(f3Prong) {
828 threeTrackArray->AddAt(negtrack1,0);
829 threeTrackArray->AddAt(postrack1,1);
830 threeTrackArray->AddAt(negtrack2,2);
dcb444c9 831 AliAODVertex* secVert3PrAOD = ReconstructSecondaryVertex(threeTrackArray,dispersion);
832 io3Prong = Make3Prong(threeTrackArray,event,secVert3PrAOD,dispersion,vertexp1n1,vertexp1n2,dcap1n1,dcap1n2,dcan1n2,ok3Prong);
b056c5e3 833 if(ok3Prong) {
2ff20727 834 AliAODVertex *v3Prong = new(verticesHFRef[iVerticesHF++])AliAODVertex(*secVert3PrAOD);
423fb9ae 835 if(!isLikeSign3Prong) {
836 rd = new(aodCharm3ProngRef[i3Prong++])AliAODRecoDecayHF3Prong(*io3Prong);
837 rd->SetSecondaryVtx(v3Prong);
838 v3Prong->SetParent(rd);
a9b75906 839 AddRefs(v3Prong,rd,event,threeTrackArray);
423fb9ae 840 } else { // isLikeSign3Prong
841 rd = new(aodLikeSign3ProngRef[iLikeSign3Prong++])AliAODRecoDecayHF3Prong(*io3Prong);
842 rd->SetSecondaryVtx(v3Prong);
843 v3Prong->SetParent(rd);
a9b75906 844 AddRefs(v3Prong,rd,event,threeTrackArray);
423fb9ae 845 }
b056c5e3 846 }
b6e66d86 847 if(io3Prong) {delete io3Prong; io3Prong=NULL;}
848 if(secVert3PrAOD) {delete secVert3PrAOD; secVert3PrAOD=NULL;}
b056c5e3 849 }
dc963de9 850
b056c5e3 851 negtrack2 = 0;
852 delete vertexp1n2;
dc963de9 853
b056c5e3 854 } // end 2nd loop on negative tracks
423fb9ae 855
b056c5e3 856 twoTrackArray2->Clear();
857
858 negtrack1 = 0;
859 delete vertexp1n1;
860 } // end 1st loop on negative tracks
861
862 postrack1 = 0;
863 } // end 1st loop on positive tracks
864
865
2ff20727 866 AliDebug(1,Form(" Total HF vertices in event = %d;",
867 (Int_t)aodVerticesHFTClArr->GetEntriesFast()));
b056c5e3 868 if(fD0toKpi) {
dcb444c9 869 AliDebug(1,Form(" D0->Kpi in event = %d;",
870 (Int_t)aodD0toKpiTClArr->GetEntriesFast()));
b056c5e3 871 }
872 if(fJPSItoEle) {
dcb444c9 873 AliDebug(1,Form(" JPSI->ee in event = %d;",
874 (Int_t)aodJPSItoEleTClArr->GetEntriesFast()));
b056c5e3 875 }
876 if(f3Prong) {
dcb444c9 877 AliDebug(1,Form(" Charm->3Prong in event = %d;",
878 (Int_t)aodCharm3ProngTClArr->GetEntriesFast()));
b056c5e3 879 }
880 if(f4Prong) {
dcb444c9 881 AliDebug(1,Form(" Charm->4Prong in event = %d;\n",
882 (Int_t)aodCharm4ProngTClArr->GetEntriesFast()));
b056c5e3 883 }
2ff20727 884 if(fDstar) {
885 AliDebug(1,Form(" D*->D0pi in event = %d;\n",
886 (Int_t)aodDstarTClArr->GetEntriesFast()));
887 }
dc963de9 888 if(fLikeSign) {
423fb9ae 889 AliDebug(1,Form(" Like-sign 2Prong in event = %d;\n",
890 (Int_t)aodLikeSign2ProngTClArr->GetEntriesFast()));
891 }
892 if(fLikeSign && f3Prong) {
893 AliDebug(1,Form(" Like-sign 3Prong in event = %d;\n",
894 (Int_t)aodLikeSign3ProngTClArr->GetEntriesFast()));
dc963de9 895 }
b056c5e3 896
897
dcb444c9 898 twoTrackArray1->Delete(); delete twoTrackArray1;
899 twoTrackArray2->Delete(); delete twoTrackArray2;
2ff20727 900 twoTrackArrayCasc->Delete(); delete twoTrackArrayCasc;
b056c5e3 901 threeTrackArray->Clear();
902 threeTrackArray->Delete(); delete threeTrackArray;
dcb444c9 903 fourTrackArray->Delete(); delete fourTrackArray;
2ff20727 904 delete [] seleFlags; seleFlags=NULL;
a25935a9 905 if(evtNumber) {delete [] evtNumber; evtNumber=NULL;}
b056c5e3 906
dcb444c9 907 if(fInputAOD) {
b6e66d86 908 seleTrksArray.Delete();
c8ab4e4f 909 if(fAODMap) { delete fAODMap; fAODMap=NULL; }
dcb444c9 910 }
b056c5e3 911
912 return;
913}
914//----------------------------------------------------------------------------
a9b75906 915void AliAnalysisVertexingHF::AddRefs(AliAODVertex *v,AliAODRecoDecayHF *rd,
916 const AliVEvent *event,
917 const TObjArray *trkArray) const
918{
919 // Add the AOD tracks as daughters of the vertex (TRef)
920 // Also add the references to the primary vertex and to the cuts
921
922 if(fInputAOD) {
923 AddDaughterRefs(v,event,trkArray);
924 rd->SetPrimaryVtxRef((AliAODVertex*)event->GetPrimaryVertex());
925 }
926
927 /*
928 rd->SetListOfCutsRef((TList*)fListOfCuts);
929 //fListOfCuts->Print();
930 cout<<fListOfCuts<<endl;
931 TList *l=(TList*)rd->GetListOfCuts();
932 cout<<l<<endl;
933 if(l) {l->Print(); }else{printf("error\n");}
934 */
935
936 return;
937}
938//----------------------------------------------------------------------------
c4cdf105 939void AliAnalysisVertexingHF::AddDaughterRefs(AliAODVertex *v,
940 const AliVEvent *event,
941 const TObjArray *trkArray) const
6a213b59 942{
dcb444c9 943 // Add the AOD tracks as daughters of the vertex (TRef)
6a213b59 944
a9b75906 945 if(v->GetNDaughters()) return; // already done
946
dcb444c9 947 Int_t nTrks = trkArray->GetEntriesFast();
6a213b59 948
e18fbfa7 949 AliExternalTrackParam *track = 0;
dcb444c9 950 AliAODTrack *aodTrack = 0;
951 Int_t id;
6a213b59 952
dcb444c9 953 for(Int_t i=0; i<nTrks; i++) {
e18fbfa7 954 track = (AliExternalTrackParam*)trkArray->UncheckedAt(i);
955 id = (Int_t)track->GetID();
2ff20727 956 //printf("---> %d\n",id);
e18fbfa7 957 if(id<0) continue; // this track is a AliAODRecoDecay
dcb444c9 958 aodTrack = (AliAODTrack*)event->GetTrack(fAODMap[id]);
959 v->AddDaughter(aodTrack);
6a213b59 960 }
6a213b59 961
962 return;
dcb444c9 963}
6a213b59 964//----------------------------------------------------------------------------
2ff20727 965AliAODRecoCascadeHF* AliAnalysisVertexingHF::MakeCascade(
966 TObjArray *twoTrackArray,AliVEvent *event,
967 AliAODVertex *secVert,
968 AliAODRecoDecayHF2Prong *rd2Prong,
969 Double_t dca,
970 Bool_t &okDstar) const
971{
972 // Make the cascade as a 2Prong decay and check if it passes Dstar
973 // reconstruction cuts
974
975 okDstar = kFALSE;
976
977 Bool_t dummy1,dummy2,dummy3;
978
979 // We use Make2Prong to construct the AliAODRecoCascadeHF
980 // (which inherits from AliAODRecoDecayHF2Prong)
981 AliAODRecoCascadeHF *theCascade =
982 (AliAODRecoCascadeHF*)Make2Prong(twoTrackArray,event,secVert,dca,
983 dummy1,dummy2,dummy3);
984 if(!theCascade) return 0x0;
985
986 // charge
987 AliESDtrack *trackPi = (AliESDtrack*)twoTrackArray->UncheckedAt(0);
988 theCascade->SetCharge(trackPi->Charge());
989
990 //--- selection cuts
991 //
992 AliAODRecoCascadeHF *tmpCascade = new AliAODRecoCascadeHF(*theCascade);
dcfa35b3 993 tmpCascade->GetSecondaryVtx()->AddDaughter(trackPi);
994 tmpCascade->GetSecondaryVtx()->AddDaughter(rd2Prong);
2ff20727 995 AliAODVertex *primVertexAOD=0;
996 if(!fRecoPrimVtxSkippingTrks && !fRmTrksFromPrimVtx) {
997 // take event primary vertex
998 primVertexAOD = PrimaryVertex();
999 tmpCascade->SetOwnPrimaryVtx(primVertexAOD);
1000 rd2Prong->SetOwnPrimaryVtx(primVertexAOD);
1001 }
1002 // select D*->D0pi
1003 if(fDstar) {
1004 Bool_t testD0=kTRUE;
1005 okDstar = tmpCascade->SelectDstar(fDstarCuts,fD0fromDstarCuts,testD0);
1006 }
dcfa35b3 1007 tmpCascade->GetSecondaryVtx()->RemoveDaughters();
2ff20727 1008 tmpCascade->UnsetOwnPrimaryVtx();
2ff20727 1009 delete tmpCascade; tmpCascade=NULL;
a25935a9 1010 if(!fRecoPrimVtxSkippingTrks && !fRmTrksFromPrimVtx && !fMixEvent) {
2ff20727 1011 rd2Prong->UnsetOwnPrimaryVtx();
2ff20727 1012 }
b6e66d86 1013 if(primVertexAOD) {delete primVertexAOD; primVertexAOD=NULL;}
2ff20727 1014 //---
1015
1016 return theCascade;
1017}
1018//-----------------------------------------------------------------------------
6a213b59 1019AliAODRecoDecayHF2Prong *AliAnalysisVertexingHF::Make2Prong(
2ff20727 1020 TObjArray *twoTrackArray,AliVEvent *event,
dcb444c9 1021 AliAODVertex *secVert,Double_t dca,
2ff20727 1022 Bool_t &okD0,Bool_t &okJPSI,
1023 Bool_t &okD0fromDstar) const
6a213b59 1024{
1025 // Make 2Prong candidates and check if they pass D0toKpi or BtoJPSI
1026 // reconstruction cuts
1027 // G.E.Bruno (J/psi), A.Dainese (D0->Kpi)
1028
2ff20727 1029 okD0=kFALSE; okJPSI=kFALSE; okD0fromDstar=kFALSE;
6a213b59 1030
1031 Double_t px[2],py[2],pz[2],d0[2],d0err[2];
2ff20727 1032 AliESDtrack *postrack = (AliESDtrack*)twoTrackArray->UncheckedAt(0);
1033 AliESDtrack *negtrack = (AliESDtrack*)twoTrackArray->UncheckedAt(1);
6a213b59 1034
1035 // propagate tracks to secondary vertex, to compute inv. mass
dcb444c9 1036 postrack->PropagateToDCA(secVert,fBzkG,kVeryBig);
1037 negtrack->PropagateToDCA(secVert,fBzkG,kVeryBig);
6a213b59 1038
1039 Double_t momentum[3];
1040 postrack->GetPxPyPz(momentum);
1041 px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2];
1042 negtrack->GetPxPyPz(momentum);
1043 px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2];
1044
1045
1046 // invariant mass cut (try to improve coding here..)
1047 Bool_t okMassCut=kFALSE;
dcb444c9 1048 if(!okMassCut && fD0toKpi) if(SelectInvMass(0,2,px,py,pz)) okMassCut=kTRUE;
6a213b59 1049 if(!okMassCut && fJPSItoEle) if(SelectInvMass(1,2,px,py,pz)) okMassCut=kTRUE;
2ff20727 1050 if(!okMassCut && fDstar) if(SelectInvMass(3,2,px,py,pz)) okMassCut=kTRUE;
6a213b59 1051 if(!okMassCut) {
dcb444c9 1052 AliDebug(2," candidate didn't pass mass cut");
6a213b59 1053 return 0x0;
1054 }
dcb444c9 1055 // primary vertex to be used by this candidate
2ff20727 1056 AliAODVertex *primVertexAOD = PrimaryVertex(twoTrackArray,event);
dcb444c9 1057 if(!primVertexAOD) return 0x0;
6a213b59 1058
a25935a9 1059
dcb444c9 1060 Double_t d0z0[2],covd0z0[3];
1061 postrack->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
6a213b59 1062 d0[0] = d0z0[0];
1063 d0err[0] = TMath::Sqrt(covd0z0[0]);
dcb444c9 1064 negtrack->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
6a213b59 1065 d0[1] = d0z0[0];
1066 d0err[1] = TMath::Sqrt(covd0z0[0]);
b6e66d86 1067
6a213b59 1068 // create the object AliAODRecoDecayHF2Prong
dcb444c9 1069 AliAODRecoDecayHF2Prong *the2Prong = new AliAODRecoDecayHF2Prong(secVert,px,py,pz,d0,d0err,dca);
6a213b59 1070 the2Prong->SetOwnPrimaryVtx(primVertexAOD);
dcb444c9 1071 UShort_t id[2]={(UShort_t)postrack->GetID(),(UShort_t)negtrack->GetID()};
1072 the2Prong->SetProngIDs(2,id);
0a65d33f 1073 delete primVertexAOD; primVertexAOD=NULL;
6a213b59 1074
f3014dc3 1075
1076 if(postrack->Charge()!=0 && negtrack->Charge()!=0) { // don't apply these cuts if it's a Dstar
1077 // select D0->Kpi
a9b75906 1078 //Int_t checkD0,checkD0bar;
1079 //if(fD0toKpi) okD0 = the2Prong->SelectD0(fD0toKpiCuts,checkD0,checkD0bar);
1080 if(fD0toKpi) okD0 = (Bool_t)fCutsD0toKpi->IsSelected(the2Prong,AliRDHFCuts::kCandidate);
f3014dc3 1081 //if(fDebug && fD0toKpi) printf(" %d\n",(Int_t)okD0);
1082 // select J/psi from B
a9b75906 1083 //Int_t checkJPSI;
1084 //if(fJPSItoEle) okJPSI = the2Prong->SelectBtoJPSI(fBtoJPSICuts,checkJPSI);
1085 if(fJPSItoEle) okJPSI = (Bool_t)fCutsJpsitoee->IsSelected(the2Prong,AliRDHFCuts::kCandidate);
f3014dc3 1086 //if(fDebug && fJPSItoEle) printf(" %d\n",(Int_t)okJPSI);
1087 // select D0->Kpi from Dstar
a9b75906 1088 //if(fDstar) okD0fromDstar = the2Prong->SelectD0(fD0fromDstarCuts,checkD0,checkD0bar);
1089 if(fDstar) okD0fromDstar = (Bool_t)fCutsD0fromDstar->IsSelected(the2Prong,AliRDHFCuts::kCandidate);
f3014dc3 1090 //if(fDebug && fDstar) printf(" %d\n",(Int_t)okD0fromDstar);
1091 }
1092
6a213b59 1093
2ff20727 1094 // remove the primary vertex (was used only for selection)
a25935a9 1095 if(!fRecoPrimVtxSkippingTrks && !fRmTrksFromPrimVtx && !fMixEvent) {
b6e66d86 1096 the2Prong->UnsetOwnPrimaryVtx();
b6e66d86 1097 }
1098
2ff20727 1099 // get PID info from ESD
5770b08b 1100 Double_t esdpid0[5]={0.,0.,0.,0.,0.};
1101 if(postrack->GetStatus()&AliESDtrack::kESDpid) postrack->GetESDpid(esdpid0);
1102 Double_t esdpid1[5]={0.,0.,0.,0.,0.};
1103 if(negtrack->GetStatus()&AliESDtrack::kESDpid) negtrack->GetESDpid(esdpid1);
2ff20727 1104 Double_t esdpid[10];
1105 for(Int_t i=0;i<5;i++) {
1106 esdpid[i] = esdpid0[i];
1107 esdpid[5+i] = esdpid1[i];
6a213b59 1108 }
2ff20727 1109 the2Prong->SetPID(2,esdpid);
6a213b59 1110
6a213b59 1111 return the2Prong;
1112}
1113//----------------------------------------------------------------------------
1114AliAODRecoDecayHF3Prong* AliAnalysisVertexingHF::Make3Prong(
dcb444c9 1115 TObjArray *threeTrackArray,AliVEvent *event,
1116 AliAODVertex *secVert,Double_t dispersion,
c4cdf105 1117 const AliAODVertex *vertexp1n1,const AliAODVertex *vertexp2n1,
6a213b59 1118 Double_t dcap1n1,Double_t dcap2n1,Double_t dcap1p2,
1119 Bool_t &ok3Prong) const
1120{
1121 // Make 3Prong candidates and check if they pass Dplus or Ds or Lambdac
1122 // reconstruction cuts
1123 // E.Bruna, F.Prino
1124
dcb444c9 1125
6a213b59 1126 ok3Prong=kFALSE;
673ac607 1127 if(!secVert || !vertexp1n1 || !vertexp2n1) return 0x0;
dcb444c9 1128
1129 Double_t px[3],py[3],pz[3],d0[3],d0err[3];
1130 Double_t momentum[3];
6a213b59 1131
6a213b59 1132
1133 AliESDtrack *postrack1 = (AliESDtrack*)threeTrackArray->UncheckedAt(0);
dcb444c9 1134 AliESDtrack *negtrack = (AliESDtrack*)threeTrackArray->UncheckedAt(1);
6a213b59 1135 AliESDtrack *postrack2 = (AliESDtrack*)threeTrackArray->UncheckedAt(2);
1136
dcb444c9 1137 postrack1->PropagateToDCA(secVert,fBzkG,kVeryBig);
1138 negtrack->PropagateToDCA(secVert,fBzkG,kVeryBig);
1139 postrack2->PropagateToDCA(secVert,fBzkG,kVeryBig);
6a213b59 1140 postrack1->GetPxPyPz(momentum);
1141 px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2];
1142 negtrack->GetPxPyPz(momentum);
1143 px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2];
1144 postrack2->GetPxPyPz(momentum);
1145 px[2] = momentum[0]; py[2] = momentum[1]; pz[2] = momentum[2];
1146
b82f6d67 1147 // invariant mass cut for D+, Ds, Lc
6a213b59 1148 Bool_t okMassCut=kFALSE;
1149 if(!okMassCut && f3Prong) if(SelectInvMass(2,3,px,py,pz)) okMassCut=kTRUE;
1150 if(!okMassCut) {
dcb444c9 1151 AliDebug(2," candidate didn't pass mass cut");
6a213b59 1152 return 0x0;
1153 }
1154
dcb444c9 1155 // primary vertex to be used by this candidate
2ff20727 1156 AliAODVertex *primVertexAOD = PrimaryVertex(threeTrackArray,event);
dcb444c9 1157 if(!primVertexAOD) return 0x0;
1158
1159 Double_t d0z0[2],covd0z0[3];
1160 postrack1->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
1161 d0[0]=d0z0[0];
1162 d0err[0] = TMath::Sqrt(covd0z0[0]);
1163 negtrack->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
1164 d0[1]=d0z0[0];
1165 d0err[1] = TMath::Sqrt(covd0z0[0]);
1166 postrack2->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
1167 d0[2]=d0z0[0];
1168 d0err[2] = TMath::Sqrt(covd0z0[0]);
6a213b59 1169
6a213b59 1170
1171 // create the object AliAODRecoDecayHF3Prong
dcb444c9 1172 Double_t pos[3]; primVertexAOD->GetXYZ(pos);
6a213b59 1173 Double_t dca[3]={dcap1n1,dcap2n1,dcap1p2};
dcb444c9 1174 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]));
1175 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]));
1176 Short_t charge=(Short_t)(postrack1->Charge()*postrack2->Charge()*negtrack->Charge());
1177 AliAODRecoDecayHF3Prong *the3Prong = new AliAODRecoDecayHF3Prong(secVert,px,py,pz,d0,d0err,dca,dispersion,dist12,dist23,charge);
6a213b59 1178 the3Prong->SetOwnPrimaryVtx(primVertexAOD);
dcb444c9 1179 UShort_t id[3]={(UShort_t)postrack1->GetID(),(UShort_t)negtrack->GetID(),(UShort_t)postrack2->GetID()};
1180 the3Prong->SetProngIDs(3,id);
6a213b59 1181
0a65d33f 1182 delete primVertexAOD; primVertexAOD=NULL;
6a213b59 1183
b82f6d67 1184 // select D+->Kpipi, Ds->KKpi, Lc->pKpi
1185 if(f3Prong) {
1186 ok3Prong = kFALSE;
a9b75906 1187 //Int_t ok1,ok2;
1188 //Int_t dum1,dum2;
1189 //if(the3Prong->SelectDplus(fDplusCuts)) ok3Prong = kTRUE;
1190 //if(the3Prong->SelectDs(fDsCuts,ok1,ok2,dum1,dum2)) ok3Prong = kTRUE;
1191 //if(the3Prong->SelectLc(fLcCuts,ok1,ok2)) ok3Prong = kTRUE;
1192
e3d40058 1193 if(fCutsDplustoKpipi->IsSelected(the3Prong,AliRDHFCuts::kCandidate)) ok3Prong = kTRUE;
1194 if(fCutsDstoKKpi->IsSelected(the3Prong,AliRDHFCuts::kCandidate)) ok3Prong = kTRUE;
1195 if(fCutsLctopKpi->IsSelected(the3Prong,AliRDHFCuts::kCandidate)) ok3Prong = kTRUE;
a9b75906 1196
b82f6d67 1197 }
6a213b59 1198 //if(fDebug) printf("ok3Prong: %d\n",(Int_t)ok3Prong);
dcb444c9 1199
a25935a9 1200 if(!fRecoPrimVtxSkippingTrks && !fRmTrksFromPrimVtx && !fMixEvent) {
b6e66d86 1201 the3Prong->UnsetOwnPrimaryVtx();
b6e66d86 1202 }
dcb444c9 1203
2ff20727 1204 // get PID info from ESD
5770b08b 1205 Double_t esdpid0[5]={0.,0.,0.,0.,0.};
1206 if(postrack1->GetStatus()&AliESDtrack::kESDpid) postrack1->GetESDpid(esdpid0);
1207 Double_t esdpid1[5]={0.,0.,0.,0.,0.};
1208 if(negtrack->GetStatus()&AliESDtrack::kESDpid) negtrack->GetESDpid(esdpid1);
1209 Double_t esdpid2[5]={0.,0.,0.,0.,0.};
1210 if(postrack2->GetStatus()&AliESDtrack::kESDpid) postrack2->GetESDpid(esdpid2);
2ff20727 1211
1212 Double_t esdpid[15];
1213 for(Int_t i=0;i<5;i++) {
1214 esdpid[i] = esdpid0[i];
1215 esdpid[5+i] = esdpid1[i];
1216 esdpid[10+i] = esdpid2[i];
6a213b59 1217 }
2ff20727 1218 the3Prong->SetPID(3,esdpid);
6a213b59 1219
6a213b59 1220 return the3Prong;
1221}
1222//----------------------------------------------------------------------------
1223AliAODRecoDecayHF4Prong* AliAnalysisVertexingHF::Make4Prong(
dcb444c9 1224 TObjArray *fourTrackArray,AliVEvent *event,
721c0b8f 1225 AliAODVertex *secVert,
c4cdf105 1226 const AliAODVertex *vertexp1n1,
1227 const AliAODVertex *vertexp1n1p2,
721c0b8f 1228 Double_t dcap1n1,Double_t dcap1n2,
1229 Double_t dcap2n1,Double_t dcap2n2,
1230 Bool_t &ok4Prong) const
6a213b59 1231{
1232 // Make 4Prong candidates and check if they pass D0toKpipipi
1233 // reconstruction cuts
1234 // G.E.Bruno, R.Romita
1235
1236 ok4Prong=kFALSE;
673ac607 1237 if(!secVert || !vertexp1n1 || !vertexp1n1p2) return 0x0;
6a213b59 1238
1239 Double_t px[4],py[4],pz[4],d0[4],d0err[4];//d0z[3];
6a213b59 1240
6a213b59 1241 AliESDtrack *postrack1 = (AliESDtrack*)fourTrackArray->UncheckedAt(0);
1242 AliESDtrack *negtrack1 = (AliESDtrack*)fourTrackArray->UncheckedAt(1);
1243 AliESDtrack *postrack2 = (AliESDtrack*)fourTrackArray->UncheckedAt(2);
1244 AliESDtrack *negtrack2 = (AliESDtrack*)fourTrackArray->UncheckedAt(3);
1245
dcb444c9 1246 postrack1->PropagateToDCA(secVert,fBzkG,kVeryBig);
1247 negtrack1->PropagateToDCA(secVert,fBzkG,kVeryBig);
1248 postrack2->PropagateToDCA(secVert,fBzkG,kVeryBig);
1249 negtrack2->PropagateToDCA(secVert,fBzkG,kVeryBig);
6a213b59 1250
1251 Double_t momentum[3];
1252 postrack1->GetPxPyPz(momentum);
1253 px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2];
1254 negtrack1->GetPxPyPz(momentum);
1255 px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2];
1256 postrack2->GetPxPyPz(momentum);
1257 px[2] = momentum[0]; py[2] = momentum[1]; pz[2] = momentum[2];
1258 negtrack2->GetPxPyPz(momentum);
1259 px[3] = momentum[0]; py[3] = momentum[1]; pz[3] = momentum[2];
1260
1261 // invariant mass cut for rho or D0 (try to improve coding here..)
0a65d33f 1262 Bool_t okMassCut=kFALSE;
c4cdf105 1263 if(!okMassCut && TMath::Abs(fD0to4ProngsCuts[8])<1.e-13){ //no PID, to be implemented with PID
0a65d33f 1264 if(SelectInvMass(4,4,px,py,pz)) okMassCut=kTRUE;
1265 }
1266 if(!okMassCut) {
1267 //if(fDebug) printf(" candidate didn't pass mass cut\n");
1268 //printf(" candidate didn't pass mass cut\n");
1269 return 0x0;
1270 }
6a213b59 1271
dcb444c9 1272 // primary vertex to be used by this candidate
2ff20727 1273 AliAODVertex *primVertexAOD = PrimaryVertex(fourTrackArray,event);
dcb444c9 1274 if(!primVertexAOD) return 0x0;
1275
721c0b8f 1276 Double_t d0z0[2],covd0z0[3];
1277 postrack1->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
1278 d0[0]=d0z0[0];
1279 d0err[0] = TMath::Sqrt(covd0z0[0]);
1280 negtrack1->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
1281 d0[1]=d0z0[0];
1282 d0err[1] = TMath::Sqrt(covd0z0[0]);
1283 postrack2->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
1284 d0[2]=d0z0[0];
1285 d0err[2] = TMath::Sqrt(covd0z0[0]);
1286 negtrack2->PropagateToDCA(primVertexAOD,fBzkG,kVeryBig,d0z0,covd0z0);
1287 d0[3]=d0z0[0];
1288 d0err[3] = TMath::Sqrt(covd0z0[0]);
1289
6a213b59 1290
1291 // create the object AliAODRecoDecayHF4Prong
dcb444c9 1292 Double_t pos[3]; primVertexAOD->GetXYZ(pos);
721c0b8f 1293 Double_t dca[6]={dcap1n1,0.,dcap1n2,dcap2n1,0.,dcap2n2};
dcb444c9 1294 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]));
721c0b8f 1295 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]));
1296 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]));
dcb444c9 1297 Short_t charge=0;
721c0b8f 1298 AliAODRecoDecayHF4Prong *the4Prong = new AliAODRecoDecayHF4Prong(secVert,px,py,pz,d0,d0err,dca,dist12,dist3,dist4,charge);
6a213b59 1299 the4Prong->SetOwnPrimaryVtx(primVertexAOD);
dcb444c9 1300 UShort_t id[4]={(UShort_t)postrack1->GetID(),(UShort_t)negtrack1->GetID(),(UShort_t)postrack2->GetID(),(UShort_t)negtrack2->GetID()};
1301 the4Prong->SetProngIDs(4,id);
6a213b59 1302
0a65d33f 1303 delete primVertexAOD; primVertexAOD=NULL;
6a213b59 1304
a9b75906 1305 //Int_t checkD0,checkD0bar;
1306 //ok4Prong=the4Prong->SelectD0(fD0to4ProngsCuts,checkD0,checkD0bar);
1307 ok4Prong=(Bool_t)fCutsD0toKpipipi->IsSelected(the4Prong,AliRDHFCuts::kCandidate);
e3d40058 1308
6a213b59 1309
a25935a9 1310 if(!fRecoPrimVtxSkippingTrks && !fRmTrksFromPrimVtx && !fMixEvent) {
b6e66d86 1311 the4Prong->UnsetOwnPrimaryVtx();
b6e66d86 1312 }
6a213b59 1313
721c0b8f 1314
6a213b59 1315 // get PID info from ESD
721c0b8f 1316 Double_t esdpid0[5]={0.,0.,0.,0.,0.};
1317 if(postrack1->GetStatus()&AliESDtrack::kESDpid) postrack1->GetESDpid(esdpid0);
1318 Double_t esdpid1[5]={0.,0.,0.,0.,0.};
1319 if(negtrack1->GetStatus()&AliESDtrack::kESDpid) negtrack1->GetESDpid(esdpid1);
1320 Double_t esdpid2[5]={0.,0.,0.,0.,0.};
1321 if(postrack2->GetStatus()&AliESDtrack::kESDpid) postrack2->GetESDpid(esdpid2);
1322 Double_t esdpid3[5]={0.,0.,0.,0.,0.};
1323 if(negtrack2->GetStatus()&AliESDtrack::kESDpid) negtrack2->GetESDpid(esdpid3);
6a213b59 1324
1325 Double_t esdpid[20];
1326 for(Int_t i=0;i<5;i++) {
2ff20727 1327 esdpid[i] = esdpid0[i];
1328 esdpid[5+i] = esdpid1[i];
6a213b59 1329 esdpid[10+i] = esdpid2[i];
1330 esdpid[15+i] = esdpid3[i];
1331 }
1332 the4Prong->SetPID(4,esdpid);
721c0b8f 1333
6a213b59 1334 return the4Prong;
1335}
1336//-----------------------------------------------------------------------------
c4cdf105 1337AliAODVertex* AliAnalysisVertexingHF::PrimaryVertex(const TObjArray *trkArray,
2ff20727 1338 AliVEvent *event) const
6a213b59 1339{
2ff20727 1340 // Returns primary vertex to be used for this candidate
6a213b59 1341
dcb444c9 1342 AliESDVertex *vertexESD = 0;
1343 AliAODVertex *vertexAOD = 0;
1344
dcb444c9 1345
1346 if(!fRecoPrimVtxSkippingTrks && !fRmTrksFromPrimVtx) {
1347 // primary vertex from the input event
1348
1349 vertexESD = new AliESDVertex(*fV1);
1350
1351 } else {
1352 // primary vertex specific to this candidate
1353
2ff20727 1354 Int_t nTrks = trkArray->GetEntriesFast();
dcb444c9 1355 AliVertexerTracks *vertexer = new AliVertexerTracks(event->GetMagneticField());
1356
1357 if(fRecoPrimVtxSkippingTrks) {
1358 // recalculating the vertex
1359
1360 if(strstr(fV1->GetTitle(),"VertexerTracksWithConstraint")) {
1361 Float_t diamondcovxy[3];
1362 event->GetDiamondCovXY(diamondcovxy);
1363 Double_t pos[3]={event->GetDiamondX(),event->GetDiamondY(),0.};
8bb09289 1364 Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
dcb444c9 1365 AliESDVertex *diamond = new AliESDVertex(pos,cov,1.,1);
1366 vertexer->SetVtxStart(diamond);
1367 delete diamond; diamond=NULL;
1368 if(strstr(fV1->GetTitle(),"VertexerTracksWithConstraintOnlyFitter"))
1369 vertexer->SetOnlyFitter();
1370 }
0a65d33f 1371 Int_t skipped[1000];
e18fbfa7 1372 Int_t nTrksToSkip=0,id;
1373 AliExternalTrackParam *t = 0;
dcb444c9 1374 for(Int_t i=0; i<nTrks; i++) {
e18fbfa7 1375 t = (AliExternalTrackParam*)trkArray->UncheckedAt(i);
1376 id = (Int_t)t->GetID();
1377 if(id<0) continue;
1378 skipped[nTrksToSkip++] = id;
dcb444c9 1379 }
0a65d33f 1380 // TEMPORARY FIX
1381 // For AOD, skip also tracks without covariance matrix
1382 if(fInputAOD) {
1383 Double_t covtest[21];
1384 for(Int_t j=0; j<event->GetNumberOfTracks(); j++) {
1385 AliVTrack *vtrack = (AliVTrack*)event->GetTrack(j);
1386 if(!vtrack->GetCovarianceXYZPxPyPz(covtest)) {
1387 id = (Int_t)t->GetID();
1388 if(id<0) continue;
1389 skipped[nTrksToSkip++] = id;
1390 }
1391 }
1392 }
1393 //
e18fbfa7 1394 vertexer->SetSkipTracks(nTrksToSkip,skipped);
dcb444c9 1395 vertexESD = (AliESDVertex*)vertexer->FindPrimaryVertex(event);
1396
1397 } else if(fRmTrksFromPrimVtx) {
1398 // removing the prongs tracks
1399
1400 TObjArray rmArray(nTrks);
1401 UShort_t *rmId = new UShort_t[nTrks];
1402 AliESDtrack *esdTrack = 0;
1403 AliESDtrack *t = 0;
1404 for(Int_t i=0; i<nTrks; i++) {
1405 t = (AliESDtrack*)trkArray->UncheckedAt(i);
1406 esdTrack = new AliESDtrack(*t);
1407 rmArray.AddLast(esdTrack);
2ff20727 1408 if(esdTrack->GetID()>=0) {
1409 rmId[i]=(UShort_t)esdTrack->GetID();
1410 } else {
1411 rmId[i]=9999;
1412 }
dcb444c9 1413 }
1414 Float_t diamondxy[2]={event->GetDiamondX(),event->GetDiamondY()};
1415 vertexESD = vertexer->RemoveTracksFromVertex(fV1,&rmArray,rmId,diamondxy);
1416 delete [] rmId; rmId=NULL;
1417 rmArray.Delete();
1418
6a213b59 1419 }
6a213b59 1420
dcb444c9 1421 if(!vertexESD) return vertexAOD;
1422 if(vertexESD->GetNContributors()<=0) {
1423 AliDebug(2,"vertexing failed");
1424 delete vertexESD; vertexESD=NULL;
1425 return vertexAOD;
6a213b59 1426 }
dcb444c9 1427
1428 delete vertexer; vertexer=NULL;
1429
6a213b59 1430 }
1431
dcb444c9 1432 // convert to AliAODVertex
1433 Double_t pos[3],cov[6],chi2perNDF;
1434 vertexESD->GetXYZ(pos); // position
1435 vertexESD->GetCovMatrix(cov); //covariance matrix
1436 chi2perNDF = vertexESD->GetChi2toNDF();
1437 delete vertexESD; vertexESD=NULL;
1438
1439 vertexAOD = new AliAODVertex(pos,cov,chi2perNDF);
6a213b59 1440
dcb444c9 1441 return vertexAOD;
6a213b59 1442}
1443//-----------------------------------------------------------------------------
1444void AliAnalysisVertexingHF::PrintStatus() const {
1445 // Print parameters being used
1446
f8fa4595 1447 //printf("Preselections:\n");
1448 // fTrackFilter->Dump();
4d5c9633 1449 if(fSecVtxWithKF) {
1450 printf("Secondary vertex with Kalman filter package (AliKFParticle)\n");
1451 } else {
1452 printf("Secondary vertex with AliVertexerTracks\n");
1453 }
6a213b59 1454 if(fRecoPrimVtxSkippingTrks) printf("RecoPrimVtxSkippingTrks\n");
1455 if(fRmTrksFromPrimVtx) printf("RmTrksFromPrimVtx\n");
1456 if(fD0toKpi) {
1457 printf("Reconstruct D0->Kpi candidates with cuts:\n");
e3d40058 1458 if(fCutsD0toKpi) fCutsD0toKpi->PrintAll();
a9b75906 1459 /*
6a213b59 1460 printf(" |M-MD0| [GeV] < %f\n",fD0toKpiCuts[0]);
1461 printf(" dca [cm] < %f\n",fD0toKpiCuts[1]);
1462 printf(" cosThetaStar < %f\n",fD0toKpiCuts[2]);
1463 printf(" pTK [GeV/c] > %f\n",fD0toKpiCuts[3]);
1464 printf(" pTpi [GeV/c] > %f\n",fD0toKpiCuts[4]);
1465 printf(" |d0K| [cm] < %f\n",fD0toKpiCuts[5]);
1466 printf(" |d0pi| [cm] < %f\n",fD0toKpiCuts[6]);
1467 printf(" d0d0 [cm^2] < %f\n",fD0toKpiCuts[7]);
1468 printf(" cosThetaPoint > %f\n",fD0toKpiCuts[8]);
a9b75906 1469 */
6a213b59 1470 }
2ff20727 1471 if(fDstar) {
eee95e18 1472 if(fFindVertexForDstar) {
1473 printf(" Reconstruct a secondary vertex for the D*\n");
1474 } else {
1475 printf(" Assume the D* comes from the primary vertex\n");
1476 }
2ff20727 1477 printf(" |M-MD*| [GeV] < %f\n",fDstarCuts[0]);
1478 printf(" |M_Kpipi-M_Kpi-(MD*-MD0)| [GeV] < %f\n",fDstarCuts[1]);
1479 printf(" pTpisoft [GeV/c] > %f\n",fDstarCuts[2]);
1480 printf(" pTpisoft [GeV/c] < %f\n",fDstarCuts[3]);
1481 printf(" Theta(pisoft,D0plane) < %f\n",fDstarCuts[4]);
1482 printf("Reconstruct D*->D0pi candidates with cuts:\n");
1483 printf(" D0 from D* cuts:\n");
e3d40058 1484 if(fCutsD0fromDstar) fCutsD0fromDstar->PrintAll();
a9b75906 1485 /*
2ff20727 1486 printf(" |M-MD0| [GeV] < %f\n",fD0fromDstarCuts[0]);
1487 printf(" dca [cm] < %f\n",fD0fromDstarCuts[1]);
1488 printf(" cosThetaStar < %f\n",fD0fromDstarCuts[2]);
1489 printf(" pTK [GeV/c] > %f\n",fD0fromDstarCuts[3]);
1490 printf(" pTpi [GeV/c] > %f\n",fD0fromDstarCuts[4]);
1491 printf(" |d0K| [cm] < %f\n",fD0fromDstarCuts[5]);
1492 printf(" |d0pi| [cm] < %f\n",fD0fromDstarCuts[6]);
1493 printf(" d0d0 [cm^2] < %f\n",fD0fromDstarCuts[7]);
1494 printf(" cosThetaPoint > %f\n",fD0fromDstarCuts[8]);
a9b75906 1495 */
2ff20727 1496 }
6a213b59 1497 if(fJPSItoEle) {
1498 printf("Reconstruct J/psi from B candidates with cuts:\n");
e3d40058 1499 if(fCutsJpsitoee) fCutsJpsitoee->PrintAll();
a9b75906 1500 /*
6a213b59 1501 printf(" |M-MJPSI| [GeV] < %f\n",fBtoJPSICuts[0]);
1502 printf(" dca [cm] < %f\n",fBtoJPSICuts[1]);
1503 printf(" cosThetaStar < %f\n",fBtoJPSICuts[2]);
1504 printf(" pTP [GeV/c] > %f\n",fBtoJPSICuts[3]);
1505 printf(" pTN [GeV/c] > %f\n",fBtoJPSICuts[4]);
1506 printf(" |d0P| [cm] < %f\n",fBtoJPSICuts[5]);
1507 printf(" |d0N| [cm] < %f\n",fBtoJPSICuts[6]);
1508 printf(" d0d0 [cm^2] < %f\n",fBtoJPSICuts[7]);
1509 printf(" cosThetaPoint > %f\n",fBtoJPSICuts[8]);
a9b75906 1510 */
6a213b59 1511 }
1512 if(f3Prong) {
1513 printf("Reconstruct 3 prong candidates.\n");
4d5c9633 1514 printf(" D+->Kpipi cuts:\n");
e3d40058 1515 if(fCutsDplustoKpipi) fCutsDplustoKpipi->PrintAll();
a9b75906 1516 /*
6a213b59 1517 printf(" |M-MD+| [GeV] < %f\n",fDplusCuts[0]);
1518 printf(" pTK [GeV/c] > %f\n",fDplusCuts[1]);
1519 printf(" pTPi [GeV/c] > %f\n",fDplusCuts[2]);
1520 printf(" |d0K| [cm] > %f\n",fDplusCuts[3]);
1521 printf(" |d0Pi| [cm] > %f\n",fDplusCuts[4]);
1522 printf(" dist12 [cm] < %f\n",fDplusCuts[5]);
1523 printf(" sigmavert [cm] < %f\n",fDplusCuts[6]);
1524 printf(" dist prim-sec [cm] > %f\n",fDplusCuts[7]);
1525 printf(" pM=Max{pT1,pT2,pT3} [GeV/c] > %f\n",fDplusCuts[8]);
1526 printf(" cosThetaPoint > %f\n",fDplusCuts[9]);
1527 printf(" Sum d0^2 [cm^2] > %f\n",fDplusCuts[10]);
1528 printf(" dca cut [cm] < %f\n",fDplusCuts[11]);
a9b75906 1529 */
4d5c9633 1530 printf(" Ds->KKpi cuts:\n");
e3d40058 1531 if(fCutsDstoKKpi) fCutsDstoKKpi->PrintAll();
a9b75906 1532 /*
4d5c9633 1533 printf(" |M-MDs| [GeV] < %f\n",fDsCuts[0]);
81679460 1534 printf(" pTK [GeV/c] > %f\n",fDsCuts[1]);
1535 printf(" pTPi [GeV/c] > %f\n",fDsCuts[2]);
1536 printf(" |d0K| [cm] > %f\n",fDsCuts[3]);
1537 printf(" |d0Pi| [cm] > %f\n",fDsCuts[4]);
1538 printf(" dist12 [cm] < %f\n",fDsCuts[5]);
1539 printf(" sigmavert [cm] < %f\n",fDsCuts[6]);
1540 printf(" dist prim-sec [cm] > %f\n",fDsCuts[7]);
1541 printf(" pM=Max{pT1,pT2,pT3} [GeV/c] > %f\n",fDsCuts[8]);
1542 printf(" cosThetaPoint > %f\n",fDsCuts[9]);
1543 printf(" Sum d0^2 [cm^2] > %f\n",fDsCuts[10]);
1544 printf(" dca cut [cm] < %f\n",fDsCuts[11]);
f7c170b4 1545 printf(" Inv. Mass phi [GeV] < %f\n",fDsCuts[12]);
1546 printf(" Inv. Mass K0* [GeV] < %f\n",fDsCuts[13]);
a9b75906 1547 */
6ea608bf 1548 printf(" Lc->pKpi cuts:\n");
e3d40058 1549 if(fCutsLctopKpi) fCutsLctopKpi->PrintAll();
a9b75906 1550 /*
6ea608bf 1551 printf(" |M-MLc| [GeV] < %f\n",fLcCuts[0]);
81679460 1552 printf(" pTP [GeV/c] > %f\n",fLcCuts[1]);
1553 printf(" pTPi and pTK [GeV/c] > %f\n",fLcCuts[2]);
1554 printf(" |d0P| [cm] > %f\n",fLcCuts[3]);
1555 printf(" |d0Pi| and |d0K| [cm] > %f\n",fLcCuts[4]);
1556 printf(" dist12 [cm] < %f\n",fLcCuts[5]);
1557 printf(" sigmavert [cm] < %f\n",fLcCuts[6]);
1558 printf(" dist prim-sec [cm] > %f\n",fLcCuts[7]);
1559 printf(" pM=Max{pT1,pT2,pT3} [GeV/c] > %f\n",fLcCuts[8]);
1560 printf(" cosThetaPoint > %f\n",fLcCuts[9]);
1561 printf(" Sum d0^2 [cm^2] > %f\n",fLcCuts[10]);
1562 printf(" dca cut [cm] < %f\n",fLcCuts[11]);
a9b75906 1563 */
6a213b59 1564 }
e3d40058 1565 if(f4Prong) {
1566 printf("Reconstruct 4 prong candidates.\n");
1567 printf(" D0->Kpipipi cuts:\n");
1568 if(fCutsD0toKpipipi) fCutsD0toKpipipi->PrintAll();
1569 }
6a213b59 1570
1571 return;
1572}
1573//-----------------------------------------------------------------------------
dcb444c9 1574AliAODVertex* AliAnalysisVertexingHF::ReconstructSecondaryVertex(TObjArray *trkArray,
dcfa35b3 1575 Double_t &dispersion,Bool_t useTRefArray) const
dcb444c9 1576{
1577 // Secondary vertex reconstruction with AliVertexerTracks or AliKFParticle
1578
1579 AliESDVertex *vertexESD = 0;
1580 AliAODVertex *vertexAOD = 0;
1581
1582 if(!fSecVtxWithKF) { // AliVertexerTracks
1583
1584 AliVertexerTracks *vertexer = new AliVertexerTracks(fBzkG);
1585 vertexer->SetVtxStart(fV1);
1586 vertexESD = (AliESDVertex*)vertexer->VertexForSelectedESDTracks(trkArray);
1587 delete vertexer; vertexer=NULL;
1588
1589 if(!vertexESD) return vertexAOD;
1590
1591 if(vertexESD->GetNContributors()!=trkArray->GetEntriesFast()) {
1592 AliDebug(2,"vertexing failed");
1593 delete vertexESD; vertexESD=NULL;
1594 return vertexAOD;
1595 }
1596
1597 } else { // Kalman Filter vertexer (AliKFParticle)
1598
1599 AliKFParticle::SetField(fBzkG);
1600
a353d2e4 1601 AliKFVertex vertexKF;
dcb444c9 1602
1603 Int_t nTrks = trkArray->GetEntriesFast();
1604 for(Int_t i=0; i<nTrks; i++) {
1605 AliESDtrack *esdTrack = (AliESDtrack*)trkArray->At(i);
1606 AliKFParticle daughterKF(*esdTrack,211);
1607 vertexKF.AddDaughter(daughterKF);
1608 }
a353d2e4 1609 vertexESD = new AliESDVertex(vertexKF.Parameters(),
1610 vertexKF.CovarianceMatrix(),
1611 vertexKF.GetChi2(),
1612 vertexKF.GetNContributors());
dcb444c9 1613
1614 }
1615
1616 // convert to AliAODVertex
1617 Double_t pos[3],cov[6],chi2perNDF;
1618 vertexESD->GetXYZ(pos); // position
1619 vertexESD->GetCovMatrix(cov); //covariance matrix
1620 chi2perNDF = vertexESD->GetChi2toNDF();
1621 dispersion = vertexESD->GetDispersion();
1622 delete vertexESD; vertexESD=NULL;
1623
dcfa35b3 1624 Int_t nprongs= (useTRefArray ? 0 : trkArray->GetEntriesFast());
1625 vertexAOD = new AliAODVertex(pos,cov,chi2perNDF,0x0,-1,AliAODVertex::kUndef,nprongs);
dcb444c9 1626
1627 return vertexAOD;
1628}
1629//-----------------------------------------------------------------------------
6a213b59 1630Bool_t AliAnalysisVertexingHF::SelectInvMass(Int_t decay,
1631 Int_t nprongs,
1632 Double_t *px,
1633 Double_t *py,
1634 Double_t *pz) const {
1635 // Check invariant mass cut
1636
1637 Short_t dummycharge=0;
1638 Double_t *dummyd0 = new Double_t[nprongs];
b056c5e3 1639 for(Int_t ip=0;ip<nprongs;ip++) dummyd0[ip]=0.;
6a213b59 1640 AliAODRecoDecay *rd = new AliAODRecoDecay(0x0,nprongs,dummycharge,px,py,pz,dummyd0);
b6e66d86 1641 delete [] dummyd0; dummyd0=NULL;
6a213b59 1642
721c0b8f 1643 UInt_t pdg2[2],pdg3[3],pdg4[4];
6a213b59 1644 Double_t mPDG,minv;
1645
1646 Bool_t retval=kFALSE;
1647 switch (decay)
1648 {
1649 case 0: // D0->Kpi
1650 pdg2[0]=211; pdg2[1]=321;
1651 mPDG=TDatabasePDG::Instance()->GetParticle(421)->Mass();
1652 minv = rd->InvMass(nprongs,pdg2);
a9b75906 1653 //if(TMath::Abs(minv-mPDG)<fD0toKpiCuts[0]) retval=kTRUE;
1654 if(TMath::Abs(minv-mPDG)<fCutsD0toKpi->GetMassCut()) retval=kTRUE;
6a213b59 1655 pdg2[0]=321; pdg2[1]=211;
1656 minv = rd->InvMass(nprongs,pdg2);
a9b75906 1657 //if(TMath::Abs(minv-mPDG)<fD0toKpiCuts[0]) retval=kTRUE;
1658 if(TMath::Abs(minv-mPDG)<fCutsD0toKpi->GetMassCut()) retval=kTRUE;
6a213b59 1659 break;
1660 case 1: // JPSI->ee
1661 pdg2[0]=11; pdg2[1]=11;
1662 mPDG=TDatabasePDG::Instance()->GetParticle(443)->Mass();
1663 minv = rd->InvMass(nprongs,pdg2);
a9b75906 1664 //if(TMath::Abs(minv-mPDG)<fBtoJPSICuts[0]) retval=kTRUE;
1665 if(TMath::Abs(minv-mPDG)<fCutsJpsitoee->GetMassCut()) retval=kTRUE;
6a213b59 1666 break;
1667 case 2: // D+->Kpipi
1668 pdg3[0]=211; pdg3[1]=321; pdg3[2]=211;
1669 mPDG=TDatabasePDG::Instance()->GetParticle(411)->Mass();
1670 minv = rd->InvMass(nprongs,pdg3);
a9b75906 1671 //if(TMath::Abs(minv-mPDG)<fDplusCuts[0]) retval=kTRUE;
1672 if(TMath::Abs(minv-mPDG)<fCutsDplustoKpipi->GetMassCut()) retval=kTRUE;
b82f6d67 1673 // Ds+->KKpi
1674 pdg3[0]=321; pdg3[1]=321; pdg3[2]=211;
1675 mPDG=TDatabasePDG::Instance()->GetParticle(431)->Mass();
1676 minv = rd->InvMass(nprongs,pdg3);
a9b75906 1677 //if(TMath::Abs(minv-mPDG)<fDsCuts[0]) retval=kTRUE;
1678 if(TMath::Abs(minv-mPDG)<fCutsDstoKKpi->GetMassCut()) retval=kTRUE;
b82f6d67 1679 pdg3[0]=211; pdg3[1]=321; pdg3[2]=321;
1680 minv = rd->InvMass(nprongs,pdg3);
a9b75906 1681 //if(TMath::Abs(minv-mPDG)<fDsCuts[0]) retval=kTRUE;
1682 if(TMath::Abs(minv-mPDG)<fCutsDstoKKpi->GetMassCut()) retval=kTRUE;
b82f6d67 1683 // Lc->pKpi
1684 pdg3[0]=2212; pdg3[1]=321; pdg3[2]=211;
1685 mPDG=TDatabasePDG::Instance()->GetParticle(4122)->Mass();
1686 minv = rd->InvMass(nprongs,pdg3);
a9b75906 1687 //if(TMath::Abs(minv-mPDG)<fLcCuts[0]) retval=kTRUE;
1688 if(TMath::Abs(minv-mPDG)<fCutsLctopKpi->GetMassCut()) retval=kTRUE;
b82f6d67 1689 pdg3[0]=211; pdg3[1]=321; pdg3[2]=2212;
1690 minv = rd->InvMass(nprongs,pdg3);
a9b75906 1691 //if(TMath::Abs(minv-mPDG)<fLcCuts[0]) retval=kTRUE;
1692 if(TMath::Abs(minv-mPDG)<fCutsLctopKpi->GetMassCut()) retval=kTRUE;
6a213b59 1693 break;
2ff20727 1694 case 3: // D*->D0pi
682639e9 1695 pdg2[0]=211; pdg2[1]=421; // in twoTrackArrayCasc we put the pion first
2ff20727 1696 mPDG=TDatabasePDG::Instance()->GetParticle(413)->Mass();
1697 minv = rd->InvMass(nprongs,pdg2);
1698 if(TMath::Abs(minv-mPDG)<fDstarCuts[0]) retval=kTRUE;
1699 break;
e3d40058 1700 case 4: // D0->Kpipipi without PID
721c0b8f 1701 pdg4[0]=321; pdg4[1]=211; pdg4[2]=211; pdg4[3]=211;
1702 mPDG=TDatabasePDG::Instance()->GetParticle(421)->Mass();
1703 minv = rd->InvMass(nprongs,pdg4);
a9b75906 1704 //if(TMath::Abs(minv-mPDG)<fD0to4ProngsCuts[0]) retval=kTRUE;
1705 if(TMath::Abs(minv-mPDG)<fCutsD0toKpipipi->GetMassCut()) retval=kTRUE;
721c0b8f 1706 pdg4[0]=211; pdg4[1]=321; pdg4[2]=211; pdg4[3]=211;
1707 mPDG=TDatabasePDG::Instance()->GetParticle(421)->Mass();
1708 minv = rd->InvMass(nprongs,pdg4);
a9b75906 1709 //if(TMath::Abs(minv-mPDG)<fD0to4ProngsCuts[0]) retval=kTRUE;
1710 if(TMath::Abs(minv-mPDG)<fCutsD0toKpipipi->GetMassCut()) retval=kTRUE;
721c0b8f 1711 pdg4[0]=211; pdg4[1]=211; pdg4[2]=321; pdg4[3]=211;
1712 mPDG=TDatabasePDG::Instance()->GetParticle(421)->Mass();
1713 minv = rd->InvMass(nprongs,pdg4);
a9b75906 1714 //if(TMath::Abs(minv-mPDG)<fD0to4ProngsCuts[0]) retval=kTRUE;
1715 if(TMath::Abs(minv-mPDG)<fCutsD0toKpipipi->GetMassCut()) retval=kTRUE;
721c0b8f 1716 pdg4[0]=211; pdg4[1]=211; pdg4[2]=211; pdg4[3]=321;
1717 mPDG=TDatabasePDG::Instance()->GetParticle(421)->Mass();
1718 minv = rd->InvMass(nprongs,pdg4);
a9b75906 1719 //if(TMath::Abs(minv-mPDG)<fD0to4ProngsCuts[0]) retval=kTRUE;
1720 if(TMath::Abs(minv-mPDG)<fCutsD0toKpipipi->GetMassCut()) retval=kTRUE;
721c0b8f 1721 break;
6a213b59 1722 default:
2ff20727 1723 printf("SelectInvMass(): wrong decay selection\n");
6a213b59 1724 break;
1725 }
1726
b6e66d86 1727 delete rd; rd=NULL;
6a213b59 1728
1729 return retval;
1730}
1731//-----------------------------------------------------------------------------
c4cdf105 1732void AliAnalysisVertexingHF::SelectTracksAndCopyVertex(const AliVEvent *event,
2ff20727 1733 TObjArray &seleTrksArray,Int_t &nSeleTrks,
a25935a9 1734 UChar_t *seleFlags,Int_t *evtNumber)
6a213b59 1735{
2ff20727 1736 // Apply single-track preselection.
1737 // Fill a TObjArray with selected tracks (for displaced vertices or
1738 // soft pion from D*). Selection flag stored in seleFlags.
dcb444c9 1739 // Create the AliESDVertex object (convert from AliAODVertex if necessary)
1740 // In case of AOD input, also fill fAODMap for track index<->ID
1741
1742 const AliVVertex *vprimary = event->GetPrimaryVertex();
1743
1744 if(fV1) { delete fV1; fV1=NULL; }
c8ab4e4f 1745 if(fAODMap) { delete fAODMap; fAODMap=NULL; }
dcb444c9 1746
1747 Int_t nindices=0;
1748 UShort_t *indices = 0;
1749 Double_t pos[3],cov[6];
1750
1751 if(!fInputAOD) { // ESD
1752 fV1 = new AliESDVertex(*((AliESDVertex*)vprimary));
1753 } else { // AOD
1754 vprimary->GetXYZ(pos);
1755 vprimary->GetCovarianceMatrix(cov);
1756 fV1 = new AliESDVertex(pos,cov,100.,100,vprimary->GetName());
1757 indices = new UShort_t[event->GetNumberOfTracks()];
c8ab4e4f 1758 fAODMapSize = 100000;
1759 fAODMap = new Int_t[fAODMapSize];
dcb444c9 1760 }
1761
1762
dcb444c9 1763 Int_t entries = (Int_t)event->GetNumberOfTracks();
2ff20727 1764 Bool_t okDisplaced=kFALSE,okSoftPi=kFALSE;
1765 nSeleTrks=0;
6a213b59 1766
dcb444c9 1767 // transfer ITS tracks from event to arrays
6a213b59 1768 for(Int_t i=0; i<entries; i++) {
a25935a9 1769 AliVTrack *track;
1770 track = (AliVTrack*)event->GetTrack(i);
dcb444c9 1771
f7c170b4 1772 // skip pure ITS SA tracks
ac325db6 1773 if(track->GetStatus()&AliESDtrack::kITSpureSA) continue;
f7c170b4 1774
65dad444 1775 // TEMPORARY: check that the cov matrix is there
0ab1bd93 1776 Double_t covtest[21];
1777 if(!track->GetCovarianceXYZPxPyPz(covtest)) continue;
1778 //
65dad444 1779
dcb444c9 1780 if(fInputAOD) {
1781 AliAODTrack *aodt = (AliAODTrack*)track;
1782 if(aodt->GetUsedForPrimVtxFit()) {
1783 indices[nindices]=aodt->GetID(); nindices++;
1784 }
1785 fAODMap[(Int_t)aodt->GetID()] = i;
1786 }
6a213b59 1787
dcb444c9 1788 AliESDtrack *esdt = 0;
b6e66d86 1789
dcb444c9 1790 if(!fInputAOD) {
1791 esdt = (AliESDtrack*)track;
1792 } else {
1793 esdt = new AliESDtrack(track);
1794 }
6a213b59 1795
1796 // single track selection
2ff20727 1797 okDisplaced=kFALSE; okSoftPi=kFALSE;
a25935a9 1798 if(fMixEvent){
1799 evtNumber[i]=((AliMixedEvent*)event)->EventIndex(i);
1800 const AliVVertex* eventVtx=((AliMixedEvent*)event)->GetEventVertex(i);
1801 Double_t vtxPos[3],primPos[3],primCov[6],trasl[3];
1802 eventVtx->GetXYZ(vtxPos);
1803 vprimary->GetXYZ(primPos);
1804 eventVtx->GetCovarianceMatrix(primCov);
1805 for(Int_t ind=0;ind<3;ind++){
1806 trasl[ind]=vtxPos[ind]-primPos[ind];
1807 }
1808
1809 Bool_t isTransl=esdt->Translate(trasl,primCov);
1810 if(!isTransl) {
1811 delete esdt;
1812 esdt = NULL;
1813 continue;
1814 }
1815 }
1816
2ff20727 1817 if(SingleTrkCuts(esdt,okDisplaced,okSoftPi)) {
1818 seleTrksArray.AddLast(esdt);
1819 seleFlags[nSeleTrks]=0;
1820 if(okDisplaced) SETBIT(seleFlags[nSeleTrks],kBitDispl);
1821 if(okSoftPi) SETBIT(seleFlags[nSeleTrks],kBitSoftPi);
1822 nSeleTrks++;
1823 } else {
dcb444c9 1824 if(fInputAOD) delete esdt;
1825 esdt = NULL;
1826 continue;
2ff20727 1827 }
6a213b59 1828
dcb444c9 1829 } // end loop on tracks
1830
1831 // primary vertex from AOD
1832 if(fInputAOD) {
1833 delete fV1; fV1=NULL;
1834 vprimary->GetXYZ(pos);
1835 vprimary->GetCovarianceMatrix(cov);
1836 Double_t chi2toNDF = vprimary->GetChi2perNDF();
1837 Int_t ncontr=nindices;
a6e0ebfe 1838 if(!strcmp(vprimary->GetTitle(),"VertexerTracksWithContraint")) ncontr += 1;
dcb444c9 1839 Double_t chi2=chi2toNDF*(2.*(Double_t)ncontr-3.);
1840 fV1 = new AliESDVertex(pos,cov,chi2,ncontr,vprimary->GetName());
1841 fV1->SetTitle(vprimary->GetTitle());
1842 fV1->SetIndices(nindices,indices);
1843 }
1844 if(indices) { delete [] indices; indices=NULL; }
1845
6a213b59 1846
1847 return;
1848}
1849//-----------------------------------------------------------------------------
1850void AliAnalysisVertexingHF::SetD0toKpiCuts(Double_t cut0,Double_t cut1,
1851 Double_t cut2,Double_t cut3,Double_t cut4,
1852 Double_t cut5,Double_t cut6,
1853 Double_t cut7,Double_t cut8)
1854{
1855 // Set the cuts for D0 selection
1856 fD0toKpiCuts[0] = cut0;
1857 fD0toKpiCuts[1] = cut1;
1858 fD0toKpiCuts[2] = cut2;
1859 fD0toKpiCuts[3] = cut3;
1860 fD0toKpiCuts[4] = cut4;
1861 fD0toKpiCuts[5] = cut5;
1862 fD0toKpiCuts[6] = cut6;
1863 fD0toKpiCuts[7] = cut7;
1864 fD0toKpiCuts[8] = cut8;
1865
1866 return;
1867}
1868//-----------------------------------------------------------------------------
1869void AliAnalysisVertexingHF::SetD0toKpiCuts(const Double_t cuts[9])
1870{
1871 // Set the cuts for D0 selection
1872
1873 for(Int_t i=0; i<9; i++) fD0toKpiCuts[i] = cuts[i];
1874
1875 return;
1876}
1877//-----------------------------------------------------------------------------
2ff20727 1878void AliAnalysisVertexingHF::SetD0fromDstarCuts(Double_t cut0,Double_t cut1,
1879 Double_t cut2,Double_t cut3,Double_t cut4,
1880 Double_t cut5,Double_t cut6,
1881 Double_t cut7,Double_t cut8)
1882{
1883 // Set the cuts for D0 from D* selection
1884 fD0fromDstarCuts[0] = cut0;
1885 fD0fromDstarCuts[1] = cut1;
1886 fD0fromDstarCuts[2] = cut2;
1887 fD0fromDstarCuts[3] = cut3;
1888 fD0fromDstarCuts[4] = cut4;
1889 fD0fromDstarCuts[5] = cut5;
1890 fD0fromDstarCuts[6] = cut6;
1891 fD0fromDstarCuts[7] = cut7;
1892 fD0fromDstarCuts[8] = cut8;
1893
1894 return;
1895}
1896//-----------------------------------------------------------------------------
1897void AliAnalysisVertexingHF::SetD0fromDstarCuts(const Double_t cuts[9])
1898{
1899 // Set the cuts for D0 from D* selection
1900
1901 for(Int_t i=0; i<9; i++) fD0fromDstarCuts[i] = cuts[i];
1902
1903 return;
1904}
1905//-----------------------------------------------------------------------------
1906void AliAnalysisVertexingHF::SetDstarCuts(Double_t cut0,Double_t cut1,
1907 Double_t cut2,Double_t cut3,
1908 Double_t cut4)
1909{
1910 // Set the cuts for D* selection
1911 fDstarCuts[0] = cut0;
1912 fDstarCuts[1] = cut1;
1913 fDstarCuts[2] = cut2;
1914 fDstarCuts[3] = cut3;
1915 fDstarCuts[4] = cut4;
1916
1917 return;
1918}
1919//-----------------------------------------------------------------------------
1920void AliAnalysisVertexingHF::SetDstarCuts(const Double_t cuts[5])
1921{
1922 // Set the cuts for D* selection
1923
1924 for(Int_t i=0; i<5; i++) fDstarCuts[i] = cuts[i];
1925
1926 return;
1927}
1928//-----------------------------------------------------------------------------
6a213b59 1929void AliAnalysisVertexingHF::SetBtoJPSICuts(Double_t cut0,Double_t cut1,
1930 Double_t cut2,Double_t cut3,Double_t cut4,
1931 Double_t cut5,Double_t cut6,
1932 Double_t cut7,Double_t cut8)
1933{
1934 // Set the cuts for J/psi from B selection
1935 fBtoJPSICuts[0] = cut0;
1936 fBtoJPSICuts[1] = cut1;
1937 fBtoJPSICuts[2] = cut2;
1938 fBtoJPSICuts[3] = cut3;
1939 fBtoJPSICuts[4] = cut4;
1940 fBtoJPSICuts[5] = cut5;
1941 fBtoJPSICuts[6] = cut6;
1942 fBtoJPSICuts[7] = cut7;
1943 fBtoJPSICuts[8] = cut8;
1944
1945 return;
1946}
1947//-----------------------------------------------------------------------------
1948void AliAnalysisVertexingHF::SetBtoJPSICuts(const Double_t cuts[9])
1949{
1950 // Set the cuts for J/psi from B selection
1951
1952 for(Int_t i=0; i<9; i++) fBtoJPSICuts[i] = cuts[i];
1953
1954 return;
1955}
1956//-----------------------------------------------------------------------------
1957void AliAnalysisVertexingHF::SetDplusCuts(Double_t cut0,Double_t cut1,
1958 Double_t cut2,Double_t cut3,Double_t cut4,
1959 Double_t cut5,Double_t cut6,
1960 Double_t cut7,Double_t cut8,
1961 Double_t cut9,Double_t cut10,Double_t cut11)
1962{
1963 // Set the cuts for Dplus->Kpipi selection
1964 fDplusCuts[0] = cut0;
1965 fDplusCuts[1] = cut1;
1966 fDplusCuts[2] = cut2;
1967 fDplusCuts[3] = cut3;
1968 fDplusCuts[4] = cut4;
1969 fDplusCuts[5] = cut5;
1970 fDplusCuts[6] = cut6;
1971 fDplusCuts[7] = cut7;
1972 fDplusCuts[8] = cut8;
1973 fDplusCuts[9] = cut9;
1974 fDplusCuts[10] = cut10;
1975 fDplusCuts[11] = cut11;
1976
1977 return;
1978}
1979//-----------------------------------------------------------------------------
1980void AliAnalysisVertexingHF::SetDplusCuts(const Double_t cuts[12])
1981{
b82f6d67 1982 // Set the cuts for Dplus->Kpipi selection
6a213b59 1983
1984 for(Int_t i=0; i<12; i++) fDplusCuts[i] = cuts[i];
1985
1986 return;
1987}
1988//-----------------------------------------------------------------------------
6ea608bf 1989void AliAnalysisVertexingHF::SetDsCuts(Double_t cut0,Double_t cut1,
f7c170b4 1990 Double_t cut2,Double_t cut3,
1991 Double_t cut4,Double_t cut5,
1992 Double_t cut6,Double_t cut7,
1993 Double_t cut8,Double_t cut9,
1994 Double_t cut10,Double_t cut11,
1995 Double_t cut12,Double_t cut13)
b82f6d67 1996{
1997 // Set the cuts for Ds->KKpi selection
1998 fDsCuts[0] = cut0;
6ea608bf 1999 fDsCuts[1] = cut1;
2000 fDsCuts[2] = cut2;
2001 fDsCuts[3] = cut3;
2002 fDsCuts[4] = cut4;
2003 fDsCuts[5] = cut5;
2004 fDsCuts[6] = cut6;
2005 fDsCuts[7] = cut7;
2006 fDsCuts[8] = cut8;
2007 fDsCuts[9] = cut9;
2008 fDsCuts[10] = cut10;
2009 fDsCuts[11] = cut11;
81679460 2010 fDsCuts[12] = cut12;
f7c170b4 2011 fDsCuts[13] = cut13;
b82f6d67 2012
2013 return;
2014}
2015//-----------------------------------------------------------------------------
81679460 2016void AliAnalysisVertexingHF::SetDsCuts(const Double_t cuts[13])
b82f6d67 2017{
2018 // Set the cuts for Ds->KKpi selection
2019
f7c170b4 2020 for(Int_t i=0; i<14; i++) fDsCuts[i] = cuts[i];
b82f6d67 2021
2022 return;
2023}
2024//-----------------------------------------------------------------------------
6ea608bf 2025void AliAnalysisVertexingHF::SetLcCuts(Double_t cut0,Double_t cut1,
2026 Double_t cut2,Double_t cut3,Double_t cut4,
2027 Double_t cut5,Double_t cut6,
2028 Double_t cut7,Double_t cut8,
2029 Double_t cut9,Double_t cut10,Double_t cut11)
b82f6d67 2030{
2031 // Set the cuts for Lc->pKpi selection
2032 fLcCuts[0] = cut0;
6ea608bf 2033 fLcCuts[1] = cut1;
2034 fLcCuts[2] = cut2;
2035 fLcCuts[3] = cut3;
2036 fLcCuts[4] = cut4;
2037 fLcCuts[5] = cut5;
2038 fLcCuts[6] = cut6;
2039 fLcCuts[7] = cut7;
2040 fLcCuts[8] = cut8;
2041 fLcCuts[9] = cut9;
2042 fLcCuts[10] = cut10;
2043 fLcCuts[11] = cut11;
b82f6d67 2044
2045 return;
2046}
2047//-----------------------------------------------------------------------------
6ea608bf 2048void AliAnalysisVertexingHF::SetLcCuts(const Double_t cuts[12])
b82f6d67 2049{
2050 // Set the cuts for Lc->pKpi selection
2051
6ea608bf 2052 for(Int_t i=0; i<12; i++) fLcCuts[i] = cuts[i];
b82f6d67 2053
2054 return;
2055}
2056//-----------------------------------------------------------------------------
721c0b8f 2057void AliAnalysisVertexingHF::SetD0to4ProngsCuts(Double_t cut0,Double_t cut1,
2058 Double_t cut2,Double_t cut3,Double_t cut4,
2059 Double_t cut5,Double_t cut6,
2060 Double_t cut7,Double_t cut8)
2061{
2062 // Set the cuts for D0->Kpipipi selection
2063
2064 fD0to4ProngsCuts[0] = cut0;
2065 fD0to4ProngsCuts[1] = cut1;
2066 fD0to4ProngsCuts[2] = cut2;
2067 fD0to4ProngsCuts[3] = cut3;
2068 fD0to4ProngsCuts[4] = cut4;
2069 fD0to4ProngsCuts[5] = cut5;
2070 fD0to4ProngsCuts[6] = cut6;
2071 fD0to4ProngsCuts[7] = cut7;
2072 fD0to4ProngsCuts[8] = cut8;
2073
2074 return;
2075}
2076//-----------------------------------------------------------------------------
2077void AliAnalysisVertexingHF::SetD0to4ProngsCuts(const Double_t cuts[9])
2078{
2079 // Set the cuts for D0->Kpipipi selection
2080
2081 for(Int_t i=0; i<9; i++) fD0to4ProngsCuts[i] = cuts[i];
2082
2083 return;
2084}
2085//-----------------------------------------------------------------------------
2ff20727 2086Bool_t AliAnalysisVertexingHF::SingleTrkCuts(AliESDtrack *trk,
2087 Bool_t &okDisplaced,Bool_t &okSoftPi) const
6a213b59 2088{
2089 // Check if track passes some kinematical cuts
6a213b59 2090
f8fa4595 2091 // this is needed to store the impact parameters
2092 trk->RelateToVertex(fV1,fBzkG,kVeryBig);
2093
2ff20727 2094 UInt_t selectInfo;
f8fa4595 2095 //
2ff20727 2096 // Track selection, displaced tracks
2097 selectInfo = 0;
f8fa4595 2098 if(fTrackFilter) {
2099 selectInfo = fTrackFilter->IsSelected(trk);
6a213b59 2100 }
2ff20727 2101 if(selectInfo) okDisplaced=kTRUE;
2102 // Track selection, soft pions
2103 selectInfo = 0;
2104 if(fDstar && fTrackFilterSoftPi) {
2105 selectInfo = fTrackFilterSoftPi->IsSelected(trk);
2106 }
2107 if(selectInfo) okSoftPi=kTRUE;
6a213b59 2108
2ff20727 2109 if(okDisplaced || okSoftPi) return kTRUE;
f8fa4595 2110
2ff20727 2111 return kFALSE;
6a213b59 2112}
2113//-----------------------------------------------------------------------------