| 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 |
| 18 | // Candidates are store as objects deriving from AliAODRecoDecay. |
| 19 | // An example of usage can be found in the macro AliVertexingHFTest.C. |
| 20 | // Can be used as a task of AliAnalysisManager by means of the interface |
| 21 | // class AliAnalysisTaskVertexingHF. |
| 22 | // |
| 23 | // Origin: E.Bruna, G.E.Bruno, A.Dainese, F.Prino, R.Romita |
| 24 | //---------------------------------------------------------------------------- |
| 25 | #include <TTree.h> |
| b82f6d67 |
26 | #include <TFile.h> |
| 6a213b59 |
27 | #include <TDatabasePDG.h> |
| 28 | #include <TString.h> |
| 5e6f195c |
29 | #include "AliESDEvent.h" |
| 6a213b59 |
30 | #include "AliVertexerTracks.h" |
| b82f6d67 |
31 | #include "AliKFParticle.h" |
| 6a213b59 |
32 | #include "AliESDVertex.h" |
| 33 | #include "AliESDv0.h" |
| 34 | #include "AliAODRecoDecay.h" |
| 35 | #include "AliAODRecoDecayHF.h" |
| 36 | #include "AliAODRecoDecayHF2Prong.h" |
| 37 | #include "AliAODRecoDecayHF3Prong.h" |
| 38 | #include "AliAODRecoDecayHF4Prong.h" |
| 39 | #include "AliAnalysisVertexingHF.h" |
| 40 | |
| 41 | ClassImp(AliAnalysisVertexingHF) |
| 42 | |
| 43 | //---------------------------------------------------------------------------- |
| 44 | AliAnalysisVertexingHF::AliAnalysisVertexingHF(): |
| b82f6d67 |
45 | fBzkG(0.), |
| 46 | fSecVtxWithKF(kFALSE), |
| 6a213b59 |
47 | fRecoPrimVtxSkippingTrks(kFALSE), |
| 48 | fRmTrksFromPrimVtx(kFALSE), |
| 49 | fV1(0x0), |
| 50 | fDebug(0), |
| b82f6d67 |
51 | fD0toKpi(kTRUE), |
| 52 | fJPSItoEle(kTRUE), |
| 53 | f3Prong(kTRUE), |
| 54 | f4Prong(kTRUE), |
| 55 | fITSrefit(kFALSE), |
| 56 | fBothSPD(kTRUE), |
| 57 | fMinITSCls(5), |
| 58 | fMinPtCut(0.), |
| 59 | fMind0rphiCut(0.) |
| 6a213b59 |
60 | { |
| 61 | // Default constructor |
| 62 | |
| 63 | SetD0toKpiCuts(); |
| 64 | SetBtoJPSICuts(); |
| 65 | SetDplusCuts(); |
| b82f6d67 |
66 | SetDsCuts(); |
| 67 | SetLcCuts(); |
| 6a213b59 |
68 | } |
| 13977a79 |
69 | //-------------------------------------------------------------------------- |
| 70 | AliAnalysisVertexingHF::AliAnalysisVertexingHF(const AliAnalysisVertexingHF &source) : |
| 71 | TNamed(source), |
| b82f6d67 |
72 | fBzkG(source.fBzkG), |
| 73 | fSecVtxWithKF(source.fSecVtxWithKF), |
| 13977a79 |
74 | fRecoPrimVtxSkippingTrks(source.fRecoPrimVtxSkippingTrks), |
| 75 | fRmTrksFromPrimVtx(source.fRmTrksFromPrimVtx), |
| 76 | fV1(source.fV1), |
| 77 | fDebug(source.fDebug), |
| 78 | fD0toKpi(source.fD0toKpi), |
| 79 | fJPSItoEle(source.fJPSItoEle), |
| 80 | f3Prong(source.f3Prong), |
| 81 | f4Prong(source.f4Prong), |
| 82 | fITSrefit(source.fITSrefit), |
| 83 | fBothSPD(source.fBothSPD), |
| 84 | fMinITSCls(source.fMinITSCls), |
| 85 | fMinPtCut(source.fMinPtCut), |
| 86 | fMind0rphiCut(source.fMind0rphiCut) |
| 87 | { |
| 88 | // |
| 89 | // Copy constructor |
| 90 | // |
| 91 | for(Int_t i=0; i<9; i++) fD0toKpiCuts[i]=source.fD0toKpiCuts[i]; |
| 92 | for(Int_t i=0; i<9; i++) fBtoJPSICuts[i]=source.fBtoJPSICuts[i]; |
| 93 | for(Int_t i=0; i<12; i++) fDplusCuts[i]=source.fDplusCuts[i]; |
| b82f6d67 |
94 | for(Int_t i=0; i<1; i++) fDsCuts[i]=source.fDsCuts[i]; |
| 95 | for(Int_t i=0; i<1; i++) fLcCuts[i]=source.fLcCuts[i]; |
| 13977a79 |
96 | } |
| 97 | //-------------------------------------------------------------------------- |
| 98 | AliAnalysisVertexingHF &AliAnalysisVertexingHF::operator=(const AliAnalysisVertexingHF &source) |
| 99 | { |
| 100 | // |
| 101 | // assignment operator |
| 102 | // |
| 103 | if(&source == this) return *this; |
| b82f6d67 |
104 | fBzkG = source.fBzkG; |
| 105 | fSecVtxWithKF = source.fSecVtxWithKF; |
| 13977a79 |
106 | fRecoPrimVtxSkippingTrks = source.fRecoPrimVtxSkippingTrks; |
| 107 | fRmTrksFromPrimVtx = source.fRmTrksFromPrimVtx; |
| 108 | fV1 = source.fV1; |
| 109 | fDebug = source.fDebug; |
| 110 | fD0toKpi = source.fD0toKpi; |
| 111 | fJPSItoEle = source.fJPSItoEle; |
| 112 | f3Prong = source.f3Prong; |
| 113 | f4Prong = source.f4Prong; |
| 114 | fITSrefit = source.fITSrefit; |
| 115 | fBothSPD = source.fBothSPD; |
| 116 | fMinITSCls = source.fMinITSCls; |
| 117 | fMinPtCut = source.fMinPtCut; |
| 118 | fMind0rphiCut = source.fMind0rphiCut; |
| 119 | |
| 120 | for(Int_t i=0; i<9; i++) fD0toKpiCuts[i]=source.fD0toKpiCuts[i]; |
| 121 | for(Int_t i=0; i<9; i++) fBtoJPSICuts[i]=source.fBtoJPSICuts[i]; |
| 122 | for(Int_t i=0; i<12; i++) fDplusCuts[i]=source.fDplusCuts[i]; |
| b82f6d67 |
123 | for(Int_t i=0; i<1; i++) fDsCuts[i]=source.fDsCuts[i]; |
| 124 | for(Int_t i=0; i<1; i++) fLcCuts[i]=source.fLcCuts[i]; |
| 13977a79 |
125 | |
| 126 | return *this; |
| 127 | } |
| 6a213b59 |
128 | //---------------------------------------------------------------------------- |
| 129 | AliAnalysisVertexingHF::~AliAnalysisVertexingHF() { |
| 130 | // Destructor |
| 13977a79 |
131 | if(fV1) { delete fV1; fV1=0; } |
| 6a213b59 |
132 | } |
| 133 | //---------------------------------------------------------------------------- |
| 260836a9 |
134 | void AliAnalysisVertexingHF::FindCandidates(AliESDEvent *esd,TTree *treeout[]) |
| 6a213b59 |
135 | { |
| 136 | // Find heavy-flavour vertex candidates |
| 6a213b59 |
137 | |
| b82f6d67 |
138 | AliAODRecoDecayHF2Prong *io2Prong = new AliAODRecoDecayHF2Prong(); |
| 139 | AliAODRecoDecayHF3Prong *io3Prong = new AliAODRecoDecayHF3Prong(); |
| 140 | AliAODRecoDecayHF4Prong *io4Prong = new AliAODRecoDecayHF4Prong(); |
| 6a213b59 |
141 | Int_t itree=0; |
| 142 | Int_t itreeD0toKpi=-1,itreeJPSItoEle=-1,itree3Prong=-1,itree4Prong=-1; |
| 143 | Int_t initEntriesD0toKpi=0,initEntriesJPSItoEle=0,initEntries3Prong=0,initEntries4Prong=0; |
| 144 | if(fD0toKpi) { |
| 145 | itreeD0toKpi=itree; |
| 260836a9 |
146 | treeout[itree]->SetBranchAddress("D0toKpi",&io2Prong); |
| 6a213b59 |
147 | itree++; |
| 260836a9 |
148 | initEntriesD0toKpi = treeout[itreeD0toKpi]->GetEntries(); |
| 6a213b59 |
149 | } |
| 150 | if(fJPSItoEle) { |
| 151 | itreeJPSItoEle=itree; |
| 260836a9 |
152 | treeout[itree]->SetBranchAddress("JPSItoEle",&io2Prong); |
| 6a213b59 |
153 | itree++; |
| 260836a9 |
154 | initEntriesJPSItoEle = treeout[itreeJPSItoEle]->GetEntries(); |
| 6a213b59 |
155 | } |
| 156 | if(f3Prong) { |
| 157 | itree3Prong=itree; |
| 260836a9 |
158 | treeout[itree]->SetBranchAddress("Charmto3Prong",&io3Prong); |
| 6a213b59 |
159 | itree++; |
| 260836a9 |
160 | initEntries3Prong = treeout[itree3Prong]->GetEntries(); |
| 6a213b59 |
161 | } |
| 162 | if(f4Prong) { |
| 163 | itree4Prong=itree; |
| 260836a9 |
164 | treeout[itree]->SetBranchAddress("D0to4Prong",&io4Prong); |
| 6a213b59 |
165 | itree++; |
| 260836a9 |
166 | initEntries4Prong = treeout[itree4Prong]->GetEntries(); |
| 6a213b59 |
167 | } |
| 168 | delete io2Prong; io2Prong = NULL; |
| 169 | delete io3Prong; io3Prong = NULL; |
| 170 | delete io4Prong; io4Prong = NULL; |
| 171 | |
| 172 | Int_t iTrkP1,iTrkP2,iTrkN1,iTrkN2,trkEntries; |
| 173 | Int_t nTrksP=0,nTrksN=0; |
| 174 | Double_t xdummy,ydummy,dcap1n1,dcap1n2,dcap2n1,dcap1p2,dcan1n2; |
| 175 | Bool_t okD0=kFALSE,okJPSI=kFALSE,ok3Prong=kFALSE,ok4Prong=kFALSE; |
| 176 | AliESDtrack *postrack1 = 0; |
| 177 | AliESDtrack *postrack2 = 0; |
| 178 | AliESDtrack *negtrack1 = 0; |
| 179 | AliESDtrack *negtrack2 = 0; |
| 180 | Double_t dcaMax = fD0toKpiCuts[1]; |
| 181 | if(dcaMax<fBtoJPSICuts[1]) dcaMax=fBtoJPSICuts[1]; |
| 182 | if(dcaMax<fDplusCuts[11]) dcaMax=fDplusCuts[11]; |
| 183 | if(fDebug) printf(" dca cut set to %f cm\n",dcaMax); |
| 184 | |
| 6a213b59 |
185 | Int_t ev = (Int_t)esd->GetEventNumberInFile(); |
| 186 | printf("--- Finding candidates in event %d\n",ev); |
| 187 | |
| b82f6d67 |
188 | fBzkG = (Double_t)esd->GetMagneticField(); |
| 6a213b59 |
189 | |
| 190 | trkEntries = (Int_t)esd->GetNumberOfTracks(); |
| 191 | printf(" Number of tracks: %d\n",trkEntries); |
| 192 | if(trkEntries<2) return; |
| 193 | |
| b82f6d67 |
194 | // retrieve primary vertex from the AliESDEvent |
| 6a213b59 |
195 | if(!esd->GetPrimaryVertex()) { |
| 196 | printf(" No vertex in AliESD\n"); |
| 197 | return; |
| 198 | } |
| 199 | AliESDVertex copy(*(esd->GetPrimaryVertex())); |
| 200 | SetPrimaryVertex(©); |
| 6a213b59 |
201 | |
| 202 | // call function which applies sigle-track selection and |
| 203 | // separetes positives and negatives |
| 204 | TObjArray trksP(trkEntries/2); // will become TClonesArray |
| 205 | TObjArray trksN(trkEntries/2); // will become TClonesArray |
| 206 | SelectTracks(esd,trksP,nTrksP, |
| 207 | trksN,nTrksN); |
| 208 | |
| 209 | printf(" Pos. tracks: %d Neg. tracks: %d\n",nTrksP,nTrksN); |
| 210 | |
| 211 | TObjArray *twoTrackArray1 = new TObjArray(2); |
| 212 | TObjArray *twoTrackArray2 = new TObjArray(2); |
| 213 | TObjArray *threeTrackArray = new TObjArray(3); |
| 214 | TObjArray *fourTrackArray = new TObjArray(4); |
| 215 | |
| 216 | // LOOP ON POSITIVE TRACKS |
| 217 | for(iTrkP1=0; iTrkP1<nTrksP; iTrkP1++) { |
| 218 | if(fDebug) if(iTrkP1%1==0) printf(" Processing positive track number %d of %d\n",iTrkP1,nTrksP); |
| 219 | // get track from track array |
| 220 | postrack1 = (AliESDtrack*)trksP.UncheckedAt(iTrkP1); |
| 221 | |
| 222 | // LOOP ON NEGATIVE TRACKS |
| 223 | for(iTrkN1=0; iTrkN1<nTrksN; iTrkN1++) { |
| 224 | if(fDebug) if(iTrkN1%1==0) printf(" Processing negative track number %d of %d\n",iTrkN1,nTrksN); |
| 225 | // get track from tracks array |
| 226 | negtrack1 = (AliESDtrack*)trksN.UncheckedAt(iTrkN1); |
| 227 | // DCA between the two tracks |
| b82f6d67 |
228 | dcap1n1 = postrack1->GetDCA(negtrack1,fBzkG,xdummy,ydummy); |
| 6a213b59 |
229 | if(dcap1n1>dcaMax) { negtrack1=0; continue; } |
| b82f6d67 |
230 | // Vertexing |
| 6a213b59 |
231 | twoTrackArray1->AddAt(postrack1,0); |
| 232 | twoTrackArray1->AddAt(negtrack1,1); |
| b82f6d67 |
233 | AliESDVertex *vertexp1n1 = ReconstructSecondaryVertex(twoTrackArray1); |
| 6a213b59 |
234 | if(vertexp1n1->GetNContributors()!=2) { |
| 235 | if(fDebug) printf("two-track vertexing failed\n"); |
| 806873d4 |
236 | twoTrackArray1->Clear(); |
| 237 | delete vertexp1n1; |
| 238 | negtrack1=0; |
| 239 | continue; |
| 6a213b59 |
240 | } |
| 241 | if(fD0toKpi || fJPSItoEle) { |
| 242 | io2Prong = Make2Prong(twoTrackArray1,esd,vertexp1n1,dcap1n1,okD0,okJPSI); |
| 260836a9 |
243 | if(okD0) treeout[itreeD0toKpi]->Fill(); |
| 244 | if(okJPSI) treeout[itreeJPSItoEle]->Fill(); |
| 6a213b59 |
245 | delete io2Prong; io2Prong=NULL; |
| 246 | } |
| 247 | |
| 248 | twoTrackArray1->Clear(); |
| 249 | if(!f3Prong && !f4Prong) { |
| 250 | negtrack1=0; |
| 251 | delete vertexp1n1; |
| 252 | continue; |
| 253 | } |
| 254 | |
| 255 | // 2nd LOOP ON POSITIVE TRACKS |
| 256 | for(iTrkP2=iTrkP1+1; iTrkP2<nTrksP; iTrkP2++) { |
| 257 | // get track from tracks array |
| 258 | postrack2 = (AliESDtrack*)trksP.UncheckedAt(iTrkP2); |
| b82f6d67 |
259 | dcap2n1 = postrack2->GetDCA(negtrack1,fBzkG,xdummy,ydummy); |
| 6a213b59 |
260 | if(dcap2n1>dcaMax) { postrack2=0; continue; } |
| b82f6d67 |
261 | dcap1p2 = postrack2->GetDCA(postrack1,fBzkG,xdummy,ydummy); |
| 6a213b59 |
262 | if(dcap1p2>dcaMax) { postrack2=0; continue; } |
| 263 | |
| b82f6d67 |
264 | // Vertexing |
| 6a213b59 |
265 | twoTrackArray2->AddAt(postrack2,0); |
| 266 | twoTrackArray2->AddAt(negtrack1,1); |
| b82f6d67 |
267 | AliESDVertex *vertexp2n1 = ReconstructSecondaryVertex(twoTrackArray2); |
| 6a213b59 |
268 | if(f3Prong) { |
| 269 | threeTrackArray->AddAt(postrack1,0); |
| 270 | threeTrackArray->AddAt(negtrack1,1); |
| 271 | threeTrackArray->AddAt(postrack2,2); |
| 272 | io3Prong = Make3Prong(threeTrackArray,esd,vertexp1n1,vertexp2n1,dcap1n1,dcap2n1,dcap1p2,ok3Prong); |
| 260836a9 |
273 | if(ok3Prong) treeout[itree3Prong]->Fill(); |
| 6a213b59 |
274 | if(io3Prong) delete io3Prong; io3Prong=NULL; |
| 275 | } |
| 276 | if(f4Prong) { |
| 277 | // 3rd LOOP ON NEGATIVE TRACKS (for 4 prong) |
| 278 | for(iTrkN2=iTrkN1+1; iTrkN2<nTrksN; iTrkN2++) { |
| 279 | // get track from tracks array |
| 280 | negtrack2 = (AliESDtrack*)trksN.UncheckedAt(iTrkN2); |
| b82f6d67 |
281 | dcap1n2 = postrack1->GetDCA(negtrack2,fBzkG,xdummy,ydummy); |
| 6a213b59 |
282 | if(dcap1n2>dcaMax) { negtrack2=0; continue; } |
| b82f6d67 |
283 | // Vertexing |
| 6a213b59 |
284 | fourTrackArray->AddAt(postrack1,0); |
| 285 | fourTrackArray->AddAt(negtrack1,1); |
| 286 | fourTrackArray->AddAt(postrack2,2); |
| 287 | fourTrackArray->AddAt(negtrack2,3); |
| 288 | io4Prong = Make4Prong(fourTrackArray,esd,vertexp1n1,vertexp2n1,dcap1n1,dcap1n2,dcap2n1,ok4Prong); |
| 260836a9 |
289 | if(ok4Prong) treeout[itree4Prong]->Fill(); |
| 6a213b59 |
290 | delete io4Prong; io4Prong=NULL; |
| 291 | fourTrackArray->Clear(); |
| 292 | negtrack2 = 0; |
| 293 | } // end loop on negative tracks |
| 294 | } |
| 295 | postrack2 = 0; |
| 296 | delete vertexp2n1; |
| 297 | } // end 2nd loop on positive tracks |
| 298 | twoTrackArray2->Clear(); |
| 299 | |
| 300 | // 2nd LOOP ON NEGATIVE TRACKS |
| 301 | for(iTrkN2=iTrkN1+1; iTrkN2<nTrksN; iTrkN2++) { |
| 302 | // get track from tracks array |
| 303 | negtrack2 = (AliESDtrack*)trksN.UncheckedAt(iTrkN2); |
| b82f6d67 |
304 | dcap1n2 = postrack1->GetDCA(negtrack2,fBzkG,xdummy,ydummy); |
| 6a213b59 |
305 | if(dcap1n2>dcaMax) { negtrack2=0; continue; } |
| b82f6d67 |
306 | dcan1n2 = negtrack1->GetDCA(negtrack2,fBzkG,xdummy,ydummy); |
| 6a213b59 |
307 | if(dcan1n2>dcaMax) { negtrack2=0; continue; } |
| 308 | |
| b82f6d67 |
309 | // Vertexing |
| 6a213b59 |
310 | twoTrackArray2->AddAt(postrack1,0); |
| 311 | twoTrackArray2->AddAt(negtrack2,1); |
| b82f6d67 |
312 | AliESDVertex *vertexp1n2 = ReconstructSecondaryVertex(twoTrackArray2); |
| 6a213b59 |
313 | if(f3Prong) { |
| 314 | threeTrackArray->AddAt(negtrack1,0); |
| 315 | threeTrackArray->AddAt(postrack1,1); |
| 316 | threeTrackArray->AddAt(negtrack2,2); |
| 317 | io3Prong = Make3Prong(threeTrackArray,esd,vertexp1n1,vertexp1n2,dcap1n1,dcap1n2,dcan1n2,ok3Prong); |
| 260836a9 |
318 | if(ok3Prong) treeout[itree3Prong]->Fill(); |
| 6a213b59 |
319 | if(io3Prong) delete io3Prong; io3Prong=NULL; |
| 320 | } |
| 321 | negtrack2 = 0; |
| 322 | delete vertexp1n2; |
| 323 | } // end 2nd loop on negative tracks |
| 324 | twoTrackArray2->Clear(); |
| 325 | |
| 326 | negtrack1 = 0; |
| 327 | delete vertexp1n1; |
| 328 | } // end 1st loop on negative tracks |
| 329 | |
| 330 | postrack1 = 0; |
| 331 | } // end 1st loop on positive tracks |
| 332 | |
| 333 | |
| 334 | |
| 6a213b59 |
335 | //printf("delete twoTr 1\n"); |
| 336 | twoTrackArray1->Delete(); delete twoTrackArray1; |
| 337 | //printf("delete twoTr 2\n"); |
| 338 | twoTrackArray2->Delete(); delete twoTrackArray2; |
| 339 | //printf("delete threeTr 1\n"); |
| 340 | threeTrackArray->Clear(); |
| 341 | threeTrackArray->Delete(); delete threeTrackArray; |
| 342 | //printf("delete fourTr 1\n"); |
| 343 | fourTrackArray->Delete(); delete fourTrackArray; |
| 344 | |
| 345 | |
| 346 | // create a copy of this class to be written to output file |
| 347 | //AliAnalysisVertexingHF *copy = (AliAnalysisVertexingHF*)this->Clone("AnalysisVertexingHF"); |
| 348 | |
| 349 | // print statistics |
| 350 | if(fD0toKpi) { |
| 351 | printf(" D0->Kpi: event %d = %d; total = %d;\n", |
| 352 | (Int_t)esd->GetEventNumberInFile(), |
| 260836a9 |
353 | (Int_t)treeout[itreeD0toKpi]->GetEntries()-initEntriesD0toKpi, |
| 354 | (Int_t)treeout[itreeD0toKpi]->GetEntries()); |
| 6a213b59 |
355 | } |
| 356 | if(fJPSItoEle) { |
| 357 | printf(" JPSI->ee: event %d = %d; total = %d;\n", |
| 358 | (Int_t)esd->GetEventNumberInFile(), |
| 260836a9 |
359 | (Int_t)treeout[itreeJPSItoEle]->GetEntries()-initEntriesJPSItoEle, |
| 360 | (Int_t)treeout[itreeJPSItoEle]->GetEntries()); |
| 6a213b59 |
361 | } |
| 362 | if(f3Prong) { |
| 363 | printf(" Charm->3Prong: event %d = %d; total = %d;\n", |
| 364 | (Int_t)esd->GetEventNumberInFile(), |
| 260836a9 |
365 | (Int_t)treeout[itree3Prong]->GetEntries()-initEntries3Prong, |
| 366 | (Int_t)treeout[itree3Prong]->GetEntries()); |
| 6a213b59 |
367 | } |
| 368 | if(f4Prong) { |
| 369 | printf(" Charm->4Prong: event %d = %d; total = %d;\n", |
| 370 | (Int_t)esd->GetEventNumberInFile(), |
| 260836a9 |
371 | (Int_t)treeout[itree4Prong]->GetEntries()-initEntries4Prong, |
| 372 | (Int_t)treeout[itree4Prong]->GetEntries()); |
| 6a213b59 |
373 | } |
| 374 | |
| 375 | |
| 376 | return; |
| 377 | } |
| 378 | //---------------------------------------------------------------------------- |
| 379 | AliAODRecoDecayHF2Prong *AliAnalysisVertexingHF::Make2Prong( |
| 5e6f195c |
380 | TObjArray *twoTrackArray1,AliESDEvent *esd, |
| 6a213b59 |
381 | AliESDVertex *secVertexESD,Double_t dca, |
| 382 | Bool_t &okD0,Bool_t &okJPSI) const |
| 383 | { |
| 384 | // Make 2Prong candidates and check if they pass D0toKpi or BtoJPSI |
| 385 | // reconstruction cuts |
| 386 | // G.E.Bruno (J/psi), A.Dainese (D0->Kpi) |
| 387 | |
| 388 | okD0=kFALSE; okJPSI=kFALSE; |
| 389 | |
| 390 | Double_t px[2],py[2],pz[2],d0[2],d0err[2]; |
| 6a213b59 |
391 | |
| 392 | AliESDtrack *postrack = (AliESDtrack*)twoTrackArray1->UncheckedAt(0); |
| 393 | AliESDtrack *negtrack = (AliESDtrack*)twoTrackArray1->UncheckedAt(1); |
| 394 | |
| 395 | // propagate tracks to secondary vertex, to compute inv. mass |
| b82f6d67 |
396 | postrack->RelateToVertex(secVertexESD,fBzkG,10.); |
| 397 | negtrack->RelateToVertex(secVertexESD,fBzkG,10.); |
| 6a213b59 |
398 | |
| 399 | Double_t momentum[3]; |
| 400 | postrack->GetPxPyPz(momentum); |
| 401 | px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2]; |
| 402 | negtrack->GetPxPyPz(momentum); |
| 403 | px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2]; |
| 404 | |
| 405 | |
| 406 | // invariant mass cut (try to improve coding here..) |
| 407 | Bool_t okMassCut=kFALSE; |
| 408 | if(!okMassCut && fD0toKpi) if(SelectInvMass(0,2,px,py,pz)) okMassCut=kTRUE; |
| 409 | if(!okMassCut && fJPSItoEle) if(SelectInvMass(1,2,px,py,pz)) okMassCut=kTRUE; |
| 410 | if(!okMassCut) { |
| 411 | if(fDebug) printf(" candidate didn't pass mass cut\n"); |
| 412 | return 0x0; |
| 413 | } |
| 414 | |
| 415 | |
| 416 | AliESDVertex *primVertex = fV1; |
| 417 | AliESDVertex *ownPrimVertex=0; |
| 418 | |
| 419 | // primary vertex from *other* tracks in the event |
| 420 | if(fRecoPrimVtxSkippingTrks || fRmTrksFromPrimVtx) { |
| 421 | ownPrimVertex = OwnPrimaryVertex(2,twoTrackArray1,esd); |
| 422 | if(!ownPrimVertex) { |
| 423 | return 0x0; |
| 424 | } else { |
| 425 | if(ownPrimVertex->GetNContributors()<2) { |
| 426 | delete ownPrimVertex; |
| 427 | return 0x0; |
| 428 | } else { |
| 429 | primVertex = ownPrimVertex; |
| 430 | } |
| 431 | } |
| 432 | } |
| 433 | |
| 434 | Float_t d0z0[2],covd0z0[3]; |
| b82f6d67 |
435 | postrack->RelateToVertex(primVertex,fBzkG,10.); |
| 6a213b59 |
436 | postrack->GetImpactParameters(d0z0,covd0z0); |
| 437 | d0[0] = d0z0[0]; |
| 438 | d0err[0] = TMath::Sqrt(covd0z0[0]); |
| b82f6d67 |
439 | negtrack->RelateToVertex(primVertex,fBzkG,10.); |
| 6a213b59 |
440 | negtrack->GetImpactParameters(d0z0,covd0z0); |
| 441 | d0[1] = d0z0[0]; |
| 442 | d0err[1] = TMath::Sqrt(covd0z0[0]); |
| 443 | |
| 444 | // create the object AliAODRecoDecayHF2Prong |
| 445 | Double_t pos[3],cov[6]; |
| 446 | secVertexESD->GetXYZ(pos); // position |
| 447 | secVertexESD->GetCovMatrix(cov); //covariance matrix |
| 448 | AliAODVertex *secVertexAOD = new AliAODVertex(pos,cov,secVertexESD->GetChi2toNDF()); |
| 449 | primVertex->GetXYZ(pos); // position |
| 450 | primVertex->GetCovMatrix(cov); //covariance matrix |
| 451 | AliAODVertex *primVertexAOD = new AliAODVertex(pos,cov,primVertex->GetChi2toNDF()); |
| 452 | AliAODRecoDecayHF2Prong *the2Prong = new AliAODRecoDecayHF2Prong(secVertexAOD,px,py,pz,d0,d0err,dca); |
| b82f6d67 |
453 | the2Prong->SetOwnSecondaryVtx(secVertexAOD);//temporary |
| 6a213b59 |
454 | the2Prong->SetOwnPrimaryVtx(primVertexAOD); |
| 455 | |
| 456 | // select D0->Kpi |
| 457 | Int_t checkD0,checkD0bar; |
| 458 | if(fD0toKpi) okD0 = the2Prong->SelectD0(fD0toKpiCuts,checkD0,checkD0bar); |
| 459 | //if(fDebug && fD0toKpi) printf(" %d\n",(Int_t)okD0); |
| 460 | // select J/psi from B |
| 461 | Int_t checkJPSI; |
| 462 | if(fJPSItoEle) okJPSI = the2Prong->SelectBtoJPSI(fBtoJPSICuts,checkJPSI); |
| 463 | //if(fDebug && fJPSItoEle) printf(" %d\n",(Int_t)okJPSI); |
| 464 | |
| 465 | |
| 466 | if(okD0 || okJPSI) { |
| 467 | // get PID info from ESD |
| 468 | Double_t esdpid0[5]; |
| 469 | postrack->GetESDpid(esdpid0); |
| 470 | Double_t esdpid1[5]; |
| 471 | negtrack->GetESDpid(esdpid1); |
| 472 | Double_t esdpid[10]; |
| 473 | for(Int_t i=0;i<5;i++) { |
| 474 | esdpid[i] = esdpid0[i]; |
| 475 | esdpid[5+i] = esdpid1[i]; |
| 476 | } |
| 477 | the2Prong->SetPID(2,esdpid); |
| 478 | } |
| 479 | |
| 480 | if(ownPrimVertex) delete ownPrimVertex; |
| b82f6d67 |
481 | |
| 6a213b59 |
482 | return the2Prong; |
| 483 | } |
| 484 | //---------------------------------------------------------------------------- |
| 485 | AliAODRecoDecayHF3Prong* AliAnalysisVertexingHF::Make3Prong( |
| 5e6f195c |
486 | TObjArray *threeTrackArray,AliESDEvent *esd, |
| 6a213b59 |
487 | AliESDVertex *vertexp1n1,AliESDVertex *vertexp2n1, |
| 488 | Double_t dcap1n1,Double_t dcap2n1,Double_t dcap1p2, |
| 489 | Bool_t &ok3Prong) const |
| 490 | { |
| 491 | // Make 3Prong candidates and check if they pass Dplus or Ds or Lambdac |
| 492 | // reconstruction cuts |
| 493 | // E.Bruna, F.Prino |
| 494 | |
| 495 | ok3Prong=kFALSE; |
| 496 | Double_t px[3],py[3],pz[3],d0[3],d0err[3];//d0z[3]; |
| 497 | Float_t d0z0[2],covd0z0[3]; |
| 498 | |
| 6a213b59 |
499 | |
| 500 | AliESDtrack *postrack1 = (AliESDtrack*)threeTrackArray->UncheckedAt(0); |
| 501 | AliESDtrack *negtrack = (AliESDtrack*)threeTrackArray->UncheckedAt(1); |
| 502 | AliESDtrack *postrack2 = (AliESDtrack*)threeTrackArray->UncheckedAt(2); |
| 503 | |
| 504 | AliESDVertex *primVertex = fV1; |
| 505 | |
| b82f6d67 |
506 | postrack1->RelateToVertex(primVertex,fBzkG,10.); |
| 507 | negtrack->RelateToVertex(primVertex,fBzkG,10.); |
| 508 | postrack2->RelateToVertex(primVertex,fBzkG,10.); |
| 6a213b59 |
509 | |
| 510 | Double_t momentum[3]; |
| 511 | postrack1->GetPxPyPz(momentum); |
| 512 | px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2]; |
| 513 | negtrack->GetPxPyPz(momentum); |
| 514 | px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2]; |
| 515 | postrack2->GetPxPyPz(momentum); |
| 516 | px[2] = momentum[0]; py[2] = momentum[1]; pz[2] = momentum[2]; |
| 517 | |
| 518 | postrack1->GetImpactParameters(d0z0,covd0z0); |
| 519 | d0[0]=d0z0[0]; |
| 520 | d0err[0] = TMath::Sqrt(covd0z0[0]); |
| 521 | negtrack->GetImpactParameters(d0z0,covd0z0); |
| 522 | d0[1]=d0z0[0]; |
| 523 | d0err[1] = TMath::Sqrt(covd0z0[0]); |
| 524 | postrack2->GetImpactParameters(d0z0,covd0z0); |
| 525 | d0[2]=d0z0[0]; |
| 526 | d0err[2] = TMath::Sqrt(covd0z0[0]); |
| 527 | |
| 528 | |
| b82f6d67 |
529 | // invariant mass cut for D+, Ds, Lc |
| 6a213b59 |
530 | Bool_t okMassCut=kFALSE; |
| 531 | if(!okMassCut && f3Prong) if(SelectInvMass(2,3,px,py,pz)) okMassCut=kTRUE; |
| 532 | if(!okMassCut) { |
| 533 | if(fDebug) printf(" candidate didn't pass mass cut\n"); |
| 534 | return 0x0; |
| 535 | } |
| 536 | |
| 537 | //charge |
| 538 | Short_t charge=(Short_t)(postrack1->GetSign()*postrack2->GetSign()*negtrack->GetSign()); |
| 539 | |
| 540 | AliESDVertex *ownPrimVertex = 0; |
| 541 | // primary vertex from *other* tracks in the event |
| 542 | if(fRecoPrimVtxSkippingTrks || fRmTrksFromPrimVtx) { |
| 543 | ownPrimVertex = OwnPrimaryVertex(3,threeTrackArray,esd); |
| 544 | if(!ownPrimVertex) { |
| 545 | return 0x0; |
| 546 | } else { |
| 547 | if(ownPrimVertex->GetNContributors()<2) { |
| 548 | delete ownPrimVertex; |
| 549 | return 0x0; |
| 550 | } else { |
| 551 | primVertex = ownPrimVertex; |
| 552 | } |
| 553 | } |
| 554 | } |
| 555 | |
| 556 | // create the object AliAODRecoDecayHF3Prong |
| b82f6d67 |
557 | AliESDVertex* secVert3Prong = ReconstructSecondaryVertex(threeTrackArray); |
| 6a213b59 |
558 | Double_t pos[3],cov[6],sigmavert; |
| 559 | secVert3Prong->GetXYZ(pos); // position |
| 560 | secVert3Prong->GetCovMatrix(cov); //covariance matrix |
| 561 | sigmavert=secVert3Prong->GetDispersion(); |
| 562 | |
| 563 | AliAODVertex *secVert3PrAOD = new AliAODVertex(pos,cov,secVert3Prong->GetChi2toNDF()); |
| 564 | primVertex->GetXYZ(pos); // position |
| 565 | primVertex->GetCovMatrix(cov); //covariance matrix |
| 566 | AliAODVertex *primVertexAOD = new AliAODVertex(pos,cov,primVertex->GetChi2toNDF()); |
| 567 | Double_t dca[3]={dcap1n1,dcap2n1,dcap1p2}; |
| 568 | |
| 569 | Double_t dist12=TMath::Sqrt((vertexp1n1->GetXv()-pos[0])*(vertexp1n1->GetXv()-pos[0])+(vertexp1n1->GetYv()-pos[1])*(vertexp1n1->GetYv()-pos[1])+(vertexp1n1->GetZv()-pos[2])*(vertexp1n1->GetZv()-pos[2])); |
| 570 | Double_t dist23=TMath::Sqrt((vertexp2n1->GetXv()-pos[0])*(vertexp2n1->GetXv()-pos[0])+(vertexp2n1->GetYv()-pos[1])*(vertexp2n1->GetYv()-pos[1])+(vertexp2n1->GetZv()-pos[2])*(vertexp2n1->GetZv()-pos[2])); |
| 571 | |
| 572 | AliAODRecoDecayHF3Prong *the3Prong = new AliAODRecoDecayHF3Prong(secVert3PrAOD,px,py,pz,d0,d0err,dca,sigmavert,dist12,dist23,charge); |
| 573 | the3Prong->SetOwnPrimaryVtx(primVertexAOD); |
| b82f6d67 |
574 | the3Prong->SetOwnSecondaryVtx(secVert3PrAOD);//temporary |
| 6a213b59 |
575 | |
| 576 | |
| b82f6d67 |
577 | // select D+->Kpipi, Ds->KKpi, Lc->pKpi |
| 578 | if(f3Prong) { |
| 579 | ok3Prong = kFALSE; |
| 580 | Int_t ok1,ok2; |
| 581 | if(the3Prong->SelectDplus(fDplusCuts)) ok3Prong = kTRUE; |
| 582 | if(the3Prong->SelectDs(fDsCuts,ok1,ok2)) ok3Prong = kTRUE; |
| 583 | if(the3Prong->SelectLc(fLcCuts,ok1,ok2)) ok3Prong = kTRUE; |
| 584 | } |
| 6a213b59 |
585 | //if(fDebug) printf("ok3Prong: %d\n",(Int_t)ok3Prong); |
| 586 | if(ok3Prong) { |
| 587 | // get PID info from ESD |
| 588 | Double_t esdpid0[5]; |
| 589 | postrack1->GetESDpid(esdpid0); |
| 590 | Double_t esdpid1[5]; |
| 591 | negtrack->GetESDpid(esdpid1); |
| 592 | Double_t esdpid2[5]; |
| 593 | postrack2->GetESDpid(esdpid2); |
| 594 | |
| 595 | |
| 596 | Double_t esdpid[15]; |
| 597 | for(Int_t i=0;i<5;i++) { |
| 598 | esdpid[i] = esdpid0[i]; |
| 599 | esdpid[5+i] = esdpid1[i]; |
| 600 | esdpid[10+i] = esdpid2[i]; |
| 601 | } |
| 602 | the3Prong->SetPID(3,esdpid); |
| 603 | } |
| 604 | |
| 605 | if(ownPrimVertex) delete ownPrimVertex; |
| 606 | |
| 607 | return the3Prong; |
| 608 | } |
| 609 | //---------------------------------------------------------------------------- |
| 610 | AliAODRecoDecayHF4Prong* AliAnalysisVertexingHF::Make4Prong( |
| 5e6f195c |
611 | TObjArray *fourTrackArray,AliESDEvent *esd, |
| 6a213b59 |
612 | AliESDVertex *vertexp1n1,AliESDVertex *vertexp2n1, |
| 613 | Double_t dcap1n1,Double_t dcap1n2,Double_t dcap2n1, |
| 614 | Bool_t &ok4Prong) const |
| 615 | { |
| 616 | // Make 4Prong candidates and check if they pass D0toKpipipi |
| 617 | // reconstruction cuts |
| 618 | // G.E.Bruno, R.Romita |
| 619 | |
| 620 | ok4Prong=kFALSE; |
| 621 | |
| 622 | Double_t px[4],py[4],pz[4],d0[4],d0err[4];//d0z[3]; |
| 623 | //Float_t d0z0[2],covd0z0[3]; |
| 624 | |
| b82f6d67 |
625 | px[0]=dcap1n1*dcap1n2*dcap2n1; // TO BE CHANGED (done just to removed compilation warning about dca... not used) |
| 6a213b59 |
626 | |
| 627 | //charge |
| 628 | Short_t charge=0; |
| 629 | |
| 630 | AliESDtrack *postrack1 = (AliESDtrack*)fourTrackArray->UncheckedAt(0); |
| 631 | AliESDtrack *negtrack1 = (AliESDtrack*)fourTrackArray->UncheckedAt(1); |
| 632 | AliESDtrack *postrack2 = (AliESDtrack*)fourTrackArray->UncheckedAt(2); |
| 633 | AliESDtrack *negtrack2 = (AliESDtrack*)fourTrackArray->UncheckedAt(3); |
| 634 | |
| 635 | AliESDVertex *primVertex = fV1; |
| 636 | |
| b82f6d67 |
637 | postrack1->RelateToVertex(primVertex,fBzkG,10.); |
| 638 | negtrack1->RelateToVertex(primVertex,fBzkG,10.); |
| 639 | postrack2->RelateToVertex(primVertex,fBzkG,10.); |
| 640 | negtrack2->RelateToVertex(primVertex,fBzkG,10.); |
| 6a213b59 |
641 | |
| 642 | Double_t momentum[3]; |
| 643 | postrack1->GetPxPyPz(momentum); |
| 644 | px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2]; |
| 645 | negtrack1->GetPxPyPz(momentum); |
| 646 | px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2]; |
| 647 | postrack2->GetPxPyPz(momentum); |
| 648 | px[2] = momentum[0]; py[2] = momentum[1]; pz[2] = momentum[2]; |
| 649 | negtrack2->GetPxPyPz(momentum); |
| 650 | px[3] = momentum[0]; py[3] = momentum[1]; pz[3] = momentum[2]; |
| 651 | |
| 652 | // invariant mass cut for rho or D0 (try to improve coding here..) |
| 653 | //Bool_t okMassCut=kFALSE; |
| 654 | //if(!okMassCut) if(SelectInvMass(2,3,px,py,pz)) okMassCut=kTRUE; |
| 655 | //if(!okMassCut) { |
| 656 | // if(fDebug) printf(" candidate didn't pass mass cut\n"); |
| 657 | // return 0x0; |
| 658 | //} |
| 659 | |
| 660 | AliESDVertex *ownPrimVertex = 0; |
| 661 | // primary vertex from *other* tracks in the event |
| 662 | if(fRecoPrimVtxSkippingTrks || fRmTrksFromPrimVtx) { |
| 663 | ownPrimVertex = OwnPrimaryVertex(4,fourTrackArray,esd); |
| 664 | if(!ownPrimVertex) { |
| 665 | return 0x0; |
| 666 | } else { |
| 667 | if(ownPrimVertex->GetNContributors()<2) { |
| 668 | delete ownPrimVertex; |
| 669 | return 0x0; |
| 670 | } else { |
| 671 | primVertex = ownPrimVertex; |
| 672 | } |
| 673 | } |
| 674 | } |
| 675 | |
| 676 | // create the object AliAODRecoDecayHF4Prong |
| b82f6d67 |
677 | AliESDVertex* secVert4Prong = ReconstructSecondaryVertex(fourTrackArray); |
| 6a213b59 |
678 | Double_t pos[3],cov[6],sigmavert; |
| 679 | secVert4Prong->GetXYZ(pos); // position |
| 680 | secVert4Prong->GetCovMatrix(cov); //covariance matrix |
| 681 | sigmavert=secVert4Prong->GetDispersion(); |
| 682 | |
| 683 | AliAODVertex *secVert4PrAOD = new AliAODVertex(pos,cov,secVert4Prong->GetChi2toNDF()); |
| 684 | primVertex->GetXYZ(pos); // position |
| 685 | primVertex->GetCovMatrix(cov); //covariance matrix |
| 686 | AliAODVertex *primVertexAOD = new AliAODVertex(pos,cov,primVertex->GetChi2toNDF()); |
| 687 | //Double_t dca[6]={dcap1n1,dcap2n1,dcap1p2,0.,0.,0.}; // |
| 688 | Double_t dca[6]={0.,0.,0.,0.,0.,0.}; // modify it |
| 689 | |
| 690 | Double_t dist12=TMath::Sqrt((vertexp1n1->GetXv()-pos[0])*(vertexp1n1->GetXv()-pos[0])+(vertexp1n1->GetYv()-pos[1])*(vertexp1n1->GetYv()-pos[1])+(vertexp1n1->GetZv()-pos[2])*(vertexp1n1->GetZv()-pos[2])); |
| 691 | Double_t dist23=TMath::Sqrt((vertexp2n1->GetXv()-pos[0])*(vertexp2n1->GetXv()-pos[0])+(vertexp2n1->GetYv()-pos[1])*(vertexp2n1->GetYv()-pos[1])+(vertexp2n1->GetZv()-pos[2])*(vertexp2n1->GetZv()-pos[2])); |
| 692 | Double_t dist14=0.; // to be implemented |
| 693 | Double_t dist34=0.; // to be implemented |
| 694 | |
| 695 | //AliAODRecoDecayHF4Prong *the4Prong = new AliAODRecoDecayHF4Prong(secVert4PrAOD,px,py,pz,d0,d0err,dca,sigmavert,dist12,dist23,charge); |
| 696 | AliAODRecoDecayHF4Prong *the4Prong = new AliAODRecoDecayHF4Prong(secVert4PrAOD,px,py,pz,d0,d0err,dca,dist12,dist23,dist14,dist34,charge); |
| 697 | the4Prong->SetOwnPrimaryVtx(primVertexAOD); |
| b82f6d67 |
698 | the4Prong->SetOwnSecondaryVtx(secVert4PrAOD);//temporary |
| 6a213b59 |
699 | |
| 700 | |
| 701 | // use the following two lines once AliAODRecoDecayHF4Prong::SelectD0 is available |
| 702 | // select D0->Kpipipi |
| 703 | //Int_t checkD0,checkD0bar; |
| 704 | // ok4Prong=the4Prong->SelectD0(fD04pCuts,checkD0,checkD0bar); |
| 705 | ok4Prong=kFALSE; //for the time being ... |
| 706 | |
| 707 | |
| 708 | // get PID info from ESD |
| 709 | Double_t esdpid0[5]; |
| 710 | postrack1->GetESDpid(esdpid0); |
| 711 | Double_t esdpid1[5]; |
| 712 | negtrack1->GetESDpid(esdpid1); |
| 713 | Double_t esdpid2[5]; |
| 714 | postrack2->GetESDpid(esdpid2); |
| 715 | Double_t esdpid3[5]; |
| 716 | negtrack2->GetESDpid(esdpid3); |
| 717 | |
| 718 | Double_t esdpid[20]; |
| 719 | for(Int_t i=0;i<5;i++) { |
| 720 | esdpid[i] = esdpid0[i]; |
| 721 | esdpid[5+i] = esdpid1[i]; |
| 722 | esdpid[10+i] = esdpid2[i]; |
| 723 | esdpid[15+i] = esdpid3[i]; |
| 724 | } |
| 725 | the4Prong->SetPID(4,esdpid); |
| 726 | |
| 727 | if(ownPrimVertex) delete ownPrimVertex; |
| 728 | |
| 729 | return the4Prong; |
| 730 | } |
| 731 | //----------------------------------------------------------------------------- |
| 732 | AliESDVertex* AliAnalysisVertexingHF::OwnPrimaryVertex(Int_t ntrks, |
| 733 | TObjArray *trkArray, |
| 5e6f195c |
734 | AliESDEvent *esd) const |
| 6a213b59 |
735 | { |
| 736 | // Returns primary vertex specific to this candidate |
| 737 | |
| 738 | AliVertexerTracks *vertexer1 = new AliVertexerTracks(esd->GetMagneticField()); |
| 739 | AliESDVertex *ownPrimVertex = 0; |
| 740 | |
| 741 | // recalculating the vertex |
| 742 | if(fRecoPrimVtxSkippingTrks) { |
| 743 | if(strstr(fV1->GetTitle(),"VertexerTracksWithConstraint")) { |
| 744 | Float_t diamondcovxy[3]; |
| 745 | esd->GetDiamondCovXY(diamondcovxy); |
| 746 | Double_t pos[3]={esd->GetDiamondX(),esd->GetDiamondY(),0.}; |
| 747 | Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.}; |
| 748 | AliESDVertex *diamond = new AliESDVertex(pos,cov,1.,1); |
| 749 | vertexer1->SetVtxStart(diamond); |
| 750 | delete diamond; diamond=NULL; |
| 751 | if(strstr(fV1->GetTitle(),"VertexerTracksWithConstraintOnlyFitter")) |
| 752 | vertexer1->SetOnlyFitter(); |
| 753 | } |
| d70e06c7 |
754 | Int_t skipped[10]; |
| 6a213b59 |
755 | AliESDtrack *t = 0; |
| 756 | for(Int_t i=0; i<ntrks; i++) { |
| 757 | t = (AliESDtrack*)trkArray->UncheckedAt(i); |
| d70e06c7 |
758 | skipped[i] = (Int_t)t->GetID(); |
| 6a213b59 |
759 | } |
| 760 | vertexer1->SetSkipTracks(ntrks,skipped); |
| 761 | ownPrimVertex = (AliESDVertex*)vertexer1->FindPrimaryVertex(esd); |
| 762 | } |
| 763 | |
| 764 | // removing the prongs tracks |
| 765 | if(fRmTrksFromPrimVtx) { |
| f65f93f3 |
766 | TObjArray rmArray(ntrks); |
| 767 | UShort_t *rmId = new UShort_t[ntrks]; |
| 6a213b59 |
768 | AliESDtrack *esdTrack = 0; |
| 6a213b59 |
769 | AliESDtrack *t = 0; |
| 770 | for(Int_t i=0; i<ntrks; i++) { |
| 771 | t = (AliESDtrack*)trkArray->UncheckedAt(i); |
| 772 | esdTrack = new AliESDtrack(*t); |
| f65f93f3 |
773 | rmArray.AddLast(esdTrack); |
| 774 | rmId[i]=(UShort_t)esdTrack->GetID(); |
| 6a213b59 |
775 | delete esdTrack; |
| 776 | } |
| 777 | Float_t diamondxy[2]={esd->GetDiamondX(),esd->GetDiamondY()}; |
| f65f93f3 |
778 | ownPrimVertex = vertexer1->RemoveTracksFromVertex(fV1,&rmArray,rmId,diamondxy); |
| 779 | delete [] rmId; rmId=NULL; |
| 780 | rmArray.Delete(); |
| 6a213b59 |
781 | } |
| 782 | |
| 783 | delete vertexer1; vertexer1=NULL; |
| 784 | |
| 785 | return ownPrimVertex; |
| 786 | } |
| 787 | //----------------------------------------------------------------------------- |
| 788 | void AliAnalysisVertexingHF::PrintStatus() const { |
| 789 | // Print parameters being used |
| 790 | |
| 791 | printf("Preselections:\n"); |
| 792 | printf(" fITSrefit = %d\n",(Int_t)fITSrefit); |
| 793 | printf(" fBothSPD = %d\n",(Int_t)fBothSPD); |
| 794 | printf(" fMinITSCls = %d\n",fMinITSCls); |
| 795 | printf(" fMinPtCut = %f GeV/c\n",fMinPtCut); |
| 796 | printf(" fMind0rphiCut = %f cm\n",fMind0rphiCut); |
| 4d5c9633 |
797 | if(fSecVtxWithKF) { |
| 798 | printf("Secondary vertex with Kalman filter package (AliKFParticle)\n"); |
| 799 | } else { |
| 800 | printf("Secondary vertex with AliVertexerTracks\n"); |
| 801 | } |
| 6a213b59 |
802 | if(fRecoPrimVtxSkippingTrks) printf("RecoPrimVtxSkippingTrks\n"); |
| 803 | if(fRmTrksFromPrimVtx) printf("RmTrksFromPrimVtx\n"); |
| 804 | if(fD0toKpi) { |
| 805 | printf("Reconstruct D0->Kpi candidates with cuts:\n"); |
| 806 | printf(" |M-MD0| [GeV] < %f\n",fD0toKpiCuts[0]); |
| 807 | printf(" dca [cm] < %f\n",fD0toKpiCuts[1]); |
| 808 | printf(" cosThetaStar < %f\n",fD0toKpiCuts[2]); |
| 809 | printf(" pTK [GeV/c] > %f\n",fD0toKpiCuts[3]); |
| 810 | printf(" pTpi [GeV/c] > %f\n",fD0toKpiCuts[4]); |
| 811 | printf(" |d0K| [cm] < %f\n",fD0toKpiCuts[5]); |
| 812 | printf(" |d0pi| [cm] < %f\n",fD0toKpiCuts[6]); |
| 813 | printf(" d0d0 [cm^2] < %f\n",fD0toKpiCuts[7]); |
| 814 | printf(" cosThetaPoint > %f\n",fD0toKpiCuts[8]); |
| 815 | } |
| 816 | if(fJPSItoEle) { |
| 817 | printf("Reconstruct J/psi from B candidates with cuts:\n"); |
| 818 | printf(" |M-MJPSI| [GeV] < %f\n",fBtoJPSICuts[0]); |
| 819 | printf(" dca [cm] < %f\n",fBtoJPSICuts[1]); |
| 820 | printf(" cosThetaStar < %f\n",fBtoJPSICuts[2]); |
| 821 | printf(" pTP [GeV/c] > %f\n",fBtoJPSICuts[3]); |
| 822 | printf(" pTN [GeV/c] > %f\n",fBtoJPSICuts[4]); |
| 823 | printf(" |d0P| [cm] < %f\n",fBtoJPSICuts[5]); |
| 824 | printf(" |d0N| [cm] < %f\n",fBtoJPSICuts[6]); |
| 825 | printf(" d0d0 [cm^2] < %f\n",fBtoJPSICuts[7]); |
| 826 | printf(" cosThetaPoint > %f\n",fBtoJPSICuts[8]); |
| 827 | } |
| 828 | if(f3Prong) { |
| 829 | printf("Reconstruct 3 prong candidates.\n"); |
| 4d5c9633 |
830 | printf(" D+->Kpipi cuts:\n"); |
| 6a213b59 |
831 | printf(" |M-MD+| [GeV] < %f\n",fDplusCuts[0]); |
| 832 | printf(" pTK [GeV/c] > %f\n",fDplusCuts[1]); |
| 833 | printf(" pTPi [GeV/c] > %f\n",fDplusCuts[2]); |
| 834 | printf(" |d0K| [cm] > %f\n",fDplusCuts[3]); |
| 835 | printf(" |d0Pi| [cm] > %f\n",fDplusCuts[4]); |
| 836 | printf(" dist12 [cm] < %f\n",fDplusCuts[5]); |
| 837 | printf(" sigmavert [cm] < %f\n",fDplusCuts[6]); |
| 838 | printf(" dist prim-sec [cm] > %f\n",fDplusCuts[7]); |
| 839 | printf(" pM=Max{pT1,pT2,pT3} [GeV/c] > %f\n",fDplusCuts[8]); |
| 840 | printf(" cosThetaPoint > %f\n",fDplusCuts[9]); |
| 841 | printf(" Sum d0^2 [cm^2] > %f\n",fDplusCuts[10]); |
| 842 | printf(" dca cut [cm] < %f\n",fDplusCuts[11]); |
| 4d5c9633 |
843 | printf(" Ds->KKpi cuts:\n"); |
| 844 | printf(" |M-MDs| [GeV] < %f\n",fDsCuts[0]); |
| 845 | printf(" Lc->pKpi cuts:\n"); |
| 846 | printf(" |M-MD+| [GeV] < %f\n",fLcCuts[0]); |
| 6a213b59 |
847 | } |
| 848 | |
| 849 | return; |
| 850 | } |
| 851 | //----------------------------------------------------------------------------- |
| 852 | Bool_t AliAnalysisVertexingHF::SelectInvMass(Int_t decay, |
| 853 | Int_t nprongs, |
| 854 | Double_t *px, |
| 855 | Double_t *py, |
| 856 | Double_t *pz) const { |
| 857 | // Check invariant mass cut |
| 858 | |
| 859 | Short_t dummycharge=0; |
| 860 | Double_t *dummyd0 = new Double_t[nprongs]; |
| 861 | AliAODRecoDecay *rd = new AliAODRecoDecay(0x0,nprongs,dummycharge,px,py,pz,dummyd0); |
| 862 | delete [] dummyd0; |
| 863 | |
| 864 | UInt_t pdg2[2],pdg3[3]; |
| 865 | Double_t mPDG,minv; |
| 866 | |
| 867 | Bool_t retval=kFALSE; |
| 868 | switch (decay) |
| 869 | { |
| 870 | case 0: // D0->Kpi |
| 871 | pdg2[0]=211; pdg2[1]=321; |
| 872 | mPDG=TDatabasePDG::Instance()->GetParticle(421)->Mass(); |
| 873 | minv = rd->InvMass(nprongs,pdg2); |
| 874 | if(TMath::Abs(minv-mPDG)<fD0toKpiCuts[0]) retval=kTRUE; |
| 875 | pdg2[0]=321; pdg2[1]=211; |
| 876 | minv = rd->InvMass(nprongs,pdg2); |
| 877 | if(TMath::Abs(minv-mPDG)<fD0toKpiCuts[0]) retval=kTRUE; |
| 878 | break; |
| 879 | case 1: // JPSI->ee |
| 880 | pdg2[0]=11; pdg2[1]=11; |
| 881 | mPDG=TDatabasePDG::Instance()->GetParticle(443)->Mass(); |
| 882 | minv = rd->InvMass(nprongs,pdg2); |
| 883 | if(TMath::Abs(minv-mPDG)<fBtoJPSICuts[0]) retval=kTRUE; |
| 884 | break; |
| 885 | case 2: // D+->Kpipi |
| 886 | pdg3[0]=211; pdg3[1]=321; pdg3[2]=211; |
| 887 | mPDG=TDatabasePDG::Instance()->GetParticle(411)->Mass(); |
| 888 | minv = rd->InvMass(nprongs,pdg3); |
| 889 | if(TMath::Abs(minv-mPDG)<fDplusCuts[0]) retval=kTRUE; |
| b82f6d67 |
890 | // Ds+->KKpi |
| 891 | pdg3[0]=321; pdg3[1]=321; pdg3[2]=211; |
| 892 | mPDG=TDatabasePDG::Instance()->GetParticle(431)->Mass(); |
| 893 | minv = rd->InvMass(nprongs,pdg3); |
| 894 | if(TMath::Abs(minv-mPDG)<fDsCuts[0]) retval=kTRUE; |
| 895 | pdg3[0]=211; pdg3[1]=321; pdg3[2]=321; |
| 896 | minv = rd->InvMass(nprongs,pdg3); |
| 897 | if(TMath::Abs(minv-mPDG)<fDsCuts[0]) retval=kTRUE; |
| 898 | // Lc->pKpi |
| 899 | pdg3[0]=2212; pdg3[1]=321; pdg3[2]=211; |
| 900 | mPDG=TDatabasePDG::Instance()->GetParticle(4122)->Mass(); |
| 901 | minv = rd->InvMass(nprongs,pdg3); |
| 902 | if(TMath::Abs(minv-mPDG)<fLcCuts[0]) retval=kTRUE; |
| 903 | pdg3[0]=211; pdg3[1]=321; pdg3[2]=2212; |
| 904 | minv = rd->InvMass(nprongs,pdg3); |
| 905 | if(TMath::Abs(minv-mPDG)<fLcCuts[0]) retval=kTRUE; |
| 6a213b59 |
906 | break; |
| 907 | default: |
| 908 | break; |
| 909 | } |
| 910 | |
| 911 | delete rd; |
| 912 | |
| 913 | return retval; |
| 914 | } |
| 915 | //----------------------------------------------------------------------------- |
| 5e6f195c |
916 | void AliAnalysisVertexingHF::SelectTracks(AliESDEvent *esd, |
| 6a213b59 |
917 | TObjArray &trksP,Int_t &nTrksP, |
| 918 | TObjArray &trksN,Int_t &nTrksN) const |
| 919 | { |
| 920 | // Fill two TObjArrays with positive and negative tracks and |
| 921 | // apply single-track preselection |
| 922 | |
| 923 | nTrksP=0,nTrksN=0; |
| 924 | |
| 925 | Int_t entries = (Int_t)esd->GetNumberOfTracks(); |
| 926 | |
| 927 | // transfer ITS tracks from ESD to arrays and to a tree |
| 928 | for(Int_t i=0; i<entries; i++) { |
| 929 | |
| 930 | AliESDtrack *esdtrack = esd->GetTrack(i); |
| 931 | UInt_t status = esdtrack->GetStatus(); |
| 932 | |
| 933 | // require refit in ITS |
| 934 | if(fITSrefit && !(status&AliESDtrack::kITSrefit)) { |
| 935 | if(fDebug) printf("track %d is not kITSrefit\n",i); |
| 936 | continue; |
| 937 | } |
| 938 | |
| 939 | // require minimum # of ITS points |
| 940 | if(esdtrack->GetNcls(0)<fMinITSCls) { |
| 941 | if(fDebug) printf("track %d has %d ITS cls\n",i,esdtrack->GetNcls(0)); |
| 942 | continue; |
| 943 | } |
| 944 | // require points on the 2 pixel layers |
| 945 | if(fBothSPD) |
| 946 | if(!TESTBIT(esdtrack->GetITSClusterMap(),0) || |
| 947 | !TESTBIT(esdtrack->GetITSClusterMap(),1)) continue; |
| 948 | |
| 949 | // single track selection |
| b82f6d67 |
950 | if(!SingleTrkCuts(*esdtrack)) continue; |
| 6a213b59 |
951 | |
| 952 | if(esdtrack->GetSign()<0) { // negative track |
| 953 | trksN.AddLast(esdtrack); |
| 954 | nTrksN++; |
| 955 | } else { // positive track |
| 956 | trksP.AddLast(esdtrack); |
| 957 | nTrksP++; |
| 958 | } |
| 959 | |
| 960 | } // loop on ESD tracks |
| 961 | |
| 962 | return; |
| 963 | } |
| 964 | //----------------------------------------------------------------------------- |
| 965 | void AliAnalysisVertexingHF::SetD0toKpiCuts(Double_t cut0,Double_t cut1, |
| 966 | Double_t cut2,Double_t cut3,Double_t cut4, |
| 967 | Double_t cut5,Double_t cut6, |
| 968 | Double_t cut7,Double_t cut8) |
| 969 | { |
| 970 | // Set the cuts for D0 selection |
| 971 | fD0toKpiCuts[0] = cut0; |
| 972 | fD0toKpiCuts[1] = cut1; |
| 973 | fD0toKpiCuts[2] = cut2; |
| 974 | fD0toKpiCuts[3] = cut3; |
| 975 | fD0toKpiCuts[4] = cut4; |
| 976 | fD0toKpiCuts[5] = cut5; |
| 977 | fD0toKpiCuts[6] = cut6; |
| 978 | fD0toKpiCuts[7] = cut7; |
| 979 | fD0toKpiCuts[8] = cut8; |
| 980 | |
| 981 | return; |
| 982 | } |
| 983 | //----------------------------------------------------------------------------- |
| 984 | void AliAnalysisVertexingHF::SetD0toKpiCuts(const Double_t cuts[9]) |
| 985 | { |
| 986 | // Set the cuts for D0 selection |
| 987 | |
| 988 | for(Int_t i=0; i<9; i++) fD0toKpiCuts[i] = cuts[i]; |
| 989 | |
| 990 | return; |
| 991 | } |
| 992 | //----------------------------------------------------------------------------- |
| 993 | void AliAnalysisVertexingHF::SetBtoJPSICuts(Double_t cut0,Double_t cut1, |
| 994 | Double_t cut2,Double_t cut3,Double_t cut4, |
| 995 | Double_t cut5,Double_t cut6, |
| 996 | Double_t cut7,Double_t cut8) |
| 997 | { |
| 998 | // Set the cuts for J/psi from B selection |
| 999 | fBtoJPSICuts[0] = cut0; |
| 1000 | fBtoJPSICuts[1] = cut1; |
| 1001 | fBtoJPSICuts[2] = cut2; |
| 1002 | fBtoJPSICuts[3] = cut3; |
| 1003 | fBtoJPSICuts[4] = cut4; |
| 1004 | fBtoJPSICuts[5] = cut5; |
| 1005 | fBtoJPSICuts[6] = cut6; |
| 1006 | fBtoJPSICuts[7] = cut7; |
| 1007 | fBtoJPSICuts[8] = cut8; |
| 1008 | |
| 1009 | return; |
| 1010 | } |
| 1011 | //----------------------------------------------------------------------------- |
| 1012 | void AliAnalysisVertexingHF::SetBtoJPSICuts(const Double_t cuts[9]) |
| 1013 | { |
| 1014 | // Set the cuts for J/psi from B selection |
| 1015 | |
| 1016 | for(Int_t i=0; i<9; i++) fBtoJPSICuts[i] = cuts[i]; |
| 1017 | |
| 1018 | return; |
| 1019 | } |
| 1020 | //----------------------------------------------------------------------------- |
| 1021 | void AliAnalysisVertexingHF::SetDplusCuts(Double_t cut0,Double_t cut1, |
| 1022 | Double_t cut2,Double_t cut3,Double_t cut4, |
| 1023 | Double_t cut5,Double_t cut6, |
| 1024 | Double_t cut7,Double_t cut8, |
| 1025 | Double_t cut9,Double_t cut10,Double_t cut11) |
| 1026 | { |
| 1027 | // Set the cuts for Dplus->Kpipi selection |
| 1028 | fDplusCuts[0] = cut0; |
| 1029 | fDplusCuts[1] = cut1; |
| 1030 | fDplusCuts[2] = cut2; |
| 1031 | fDplusCuts[3] = cut3; |
| 1032 | fDplusCuts[4] = cut4; |
| 1033 | fDplusCuts[5] = cut5; |
| 1034 | fDplusCuts[6] = cut6; |
| 1035 | fDplusCuts[7] = cut7; |
| 1036 | fDplusCuts[8] = cut8; |
| 1037 | fDplusCuts[9] = cut9; |
| 1038 | fDplusCuts[10] = cut10; |
| 1039 | fDplusCuts[11] = cut11; |
| 1040 | |
| 1041 | return; |
| 1042 | } |
| 1043 | //----------------------------------------------------------------------------- |
| 1044 | void AliAnalysisVertexingHF::SetDplusCuts(const Double_t cuts[12]) |
| 1045 | { |
| b82f6d67 |
1046 | // Set the cuts for Dplus->Kpipi selection |
| 6a213b59 |
1047 | |
| 1048 | for(Int_t i=0; i<12; i++) fDplusCuts[i] = cuts[i]; |
| 1049 | |
| 1050 | return; |
| 1051 | } |
| 1052 | //----------------------------------------------------------------------------- |
| b82f6d67 |
1053 | void AliAnalysisVertexingHF::SetDsCuts(Double_t cut0) |
| 1054 | { |
| 1055 | // Set the cuts for Ds->KKpi selection |
| 1056 | fDsCuts[0] = cut0; |
| 1057 | |
| 1058 | return; |
| 1059 | } |
| 1060 | //----------------------------------------------------------------------------- |
| 1061 | void AliAnalysisVertexingHF::SetDsCuts(const Double_t cuts[1]) |
| 1062 | { |
| 1063 | // Set the cuts for Ds->KKpi selection |
| 1064 | |
| 1065 | for(Int_t i=0; i<1; i++) fDsCuts[i] = cuts[i]; |
| 1066 | |
| 1067 | return; |
| 1068 | } |
| 1069 | //----------------------------------------------------------------------------- |
| 1070 | void AliAnalysisVertexingHF::SetLcCuts(Double_t cut0) |
| 1071 | { |
| 1072 | // Set the cuts for Lc->pKpi selection |
| 1073 | fLcCuts[0] = cut0; |
| 1074 | |
| 1075 | return; |
| 1076 | } |
| 1077 | //----------------------------------------------------------------------------- |
| 1078 | void AliAnalysisVertexingHF::SetLcCuts(const Double_t cuts[1]) |
| 1079 | { |
| 1080 | // Set the cuts for Lc->pKpi selection |
| 1081 | |
| 1082 | for(Int_t i=0; i<1; i++) fLcCuts[i] = cuts[i]; |
| 1083 | |
| 1084 | return; |
| 1085 | } |
| 1086 | //----------------------------------------------------------------------------- |
| 1087 | Bool_t AliAnalysisVertexingHF::SingleTrkCuts(AliESDtrack& trk) const |
| 6a213b59 |
1088 | { |
| 1089 | // Check if track passes some kinematical cuts |
| 6a213b59 |
1090 | |
| 1091 | if(TMath::Abs(1./trk.GetParameter()[4]) < fMinPtCut) { |
| 1092 | printf("pt %f\n",1./trk.GetParameter()[4]); |
| 1093 | return kFALSE; |
| 1094 | } |
| b82f6d67 |
1095 | trk.RelateToVertex(fV1,fBzkG,10.); |
| 6a213b59 |
1096 | Float_t d0z0[2],covd0z0[3]; |
| 1097 | trk.GetImpactParameters(d0z0,covd0z0); |
| 1098 | if(TMath::Abs(d0z0[0]) < fMind0rphiCut) { |
| 1099 | printf("d0rphi %f\n",TMath::Abs(d0z0[0])); |
| 1100 | return kFALSE; |
| 1101 | } |
| 1102 | |
| 1103 | return kTRUE; |
| 1104 | } |
| 1105 | //----------------------------------------------------------------------------- |
| b82f6d67 |
1106 | AliESDVertex* AliAnalysisVertexingHF::ReconstructSecondaryVertex(TObjArray *trkArray) const |
| 1107 | { |
| 1108 | // Secondary vertex reconstruction with AliVertexerTracks or AliKFParticle |
| 1109 | |
| 1110 | AliESDVertex *vertex = 0; |
| 1111 | |
| 1112 | if(!fSecVtxWithKF) { // AliVertexerTracks |
| 1113 | |
| 1114 | AliVertexerTracks *vertexer2 = new AliVertexerTracks(fBzkG); |
| 1115 | vertexer2->SetDebug(0); |
| 1116 | vertexer2->SetVtxStart(fV1); |
| f65f93f3 |
1117 | vertex = (AliESDVertex*)vertexer2->VertexForSelectedESDTracks(trkArray); |
| b82f6d67 |
1118 | delete vertexer2; |
| 1119 | |
| 1120 | } else { // Kalman Filter vertexer (AliKFParticle) |
| 1121 | |
| 1122 | AliKFParticle::SetField(fBzkG); |
| 1123 | |
| 1124 | AliKFParticle vertexKF; |
| 1125 | |
| 1126 | Int_t nTrks = trkArray->GetEntriesFast(); |
| 1127 | for(Int_t i=0; i<nTrks; i++) { |
| 1128 | AliESDtrack *esdTrack = (AliESDtrack*)trkArray->At(i); |
| 1129 | AliKFParticle daughterKF(*esdTrack,211); |
| 1130 | vertexKF.AddDaughter(daughterKF); |
| 1131 | } |
| 1132 | vertex = new AliESDVertex(); |
| 1133 | vertexKF.CopyToESDVertex(*vertex); |
| 1134 | |
| 1135 | } |
| 1136 | |
| 1137 | return vertex; |
| 1138 | } |
| 1139 | //----------------------------------------------------------------------------- |
| 6a213b59 |
1140 | |
| 1141 | |
| 1142 | |
git://git.uio.no / u / mrichter / AliRoot.git / blame