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