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