]>
Commit | Line | Data |
---|---|---|
3244eeed | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-2006, 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 | ||
27de2dfb | 16 | /* $Id$ */ |
17 | ||
3244eeed | 18 | ///////////////////////////////////////////////////////////// |
19 | // | |
20 | // Base class for AOD reconstructed heavy-flavour decay | |
21 | // | |
22 | // Author: A.Dainese, andrea.dainese@lnl.infn.it | |
23 | ///////////////////////////////////////////////////////////// | |
24 | ||
25 | #include <TDatabasePDG.h> | |
26 | #include <TVector3.h> | |
15c6df42 | 27 | #include <TRandom.h> |
3244eeed | 28 | #include "AliAODRecoDecay.h" |
29 | #include "AliAODRecoDecayHF.h" | |
ead561f0 | 30 | #include "AliAODEvent.h" |
31 | #include "AliVertexerTracks.h" | |
32 | #include "AliExternalTrackParam.h" | |
ec653946 | 33 | #include "AliKFVertex.h" |
34 | #include "AliVVertex.h" | |
35 | #include "AliESDVertex.h" | |
3244eeed | 36 | |
37 | ClassImp(AliAODRecoDecayHF) | |
38 | ||
39 | //-------------------------------------------------------------------------- | |
40 | AliAODRecoDecayHF::AliAODRecoDecayHF() : | |
41 | AliAODRecoDecay(), | |
42 | fOwnPrimaryVtx(0x0), | |
460cd990 | 43 | fEventPrimaryVtx(), |
a9b75906 | 44 | fListOfCuts(), |
6185d025 | 45 | fd0err(0x0), |
939850df | 46 | fProngID(0x0), |
47 | fSelectionMap(0) | |
3244eeed | 48 | { |
49 | // | |
50 | // Default Constructor | |
51 | // | |
52 | } | |
53 | //-------------------------------------------------------------------------- | |
54 | AliAODRecoDecayHF::AliAODRecoDecayHF(AliAODVertex *vtx2,Int_t nprongs,Short_t charge, | |
55 | Double_t *px,Double_t *py,Double_t *pz, | |
56 | Double_t *d0,Double_t *d0err) : | |
57 | AliAODRecoDecay(vtx2,nprongs,charge,px,py,pz,d0), | |
58 | fOwnPrimaryVtx(0x0), | |
460cd990 | 59 | fEventPrimaryVtx(), |
a9b75906 | 60 | fListOfCuts(), |
6185d025 | 61 | fd0err(0x0), |
939850df | 62 | fProngID(0x0), |
63 | fSelectionMap(0) | |
3244eeed | 64 | { |
65 | // | |
66 | // Constructor with AliAODVertex for decay vertex | |
67 | // | |
68 | fd0err = new Double_t[GetNProngs()]; | |
69 | for(Int_t i=0; i<GetNProngs(); i++) fd0err[i] = d0err[i]; | |
70 | } | |
71 | //-------------------------------------------------------------------------- | |
72 | AliAODRecoDecayHF::AliAODRecoDecayHF(AliAODVertex *vtx2,Int_t nprongs,Short_t charge, | |
73 | Double_t *d0,Double_t *d0err) : | |
74 | AliAODRecoDecay(vtx2,nprongs,charge,d0), | |
75 | fOwnPrimaryVtx(0x0), | |
460cd990 | 76 | fEventPrimaryVtx(), |
a9b75906 | 77 | fListOfCuts(), |
6185d025 | 78 | fd0err(0x0), |
939850df | 79 | fProngID(0x0), |
80 | fSelectionMap(0) | |
3244eeed | 81 | { |
82 | // | |
83 | // Constructor with AliAODVertex for decay vertex and without prongs momenta | |
84 | // | |
85 | fd0err = new Double_t[GetNProngs()]; | |
86 | for(Int_t i=0; i<GetNProngs(); i++) fd0err[i] = d0err[i]; | |
87 | } | |
88 | //-------------------------------------------------------------------------- | |
b39168f9 | 89 | AliAODRecoDecayHF::AliAODRecoDecayHF(Double_t vtx1[3],Double_t vtx2[3], |
90 | Int_t nprongs,Short_t charge, | |
91 | Double_t *px,Double_t *py,Double_t *pz, | |
92 | Double_t *d0) : | |
93 | AliAODRecoDecay(0x0,nprongs,charge,px,py,pz,d0), | |
94 | fOwnPrimaryVtx(0x0), | |
460cd990 | 95 | fEventPrimaryVtx(), |
a9b75906 | 96 | fListOfCuts(), |
b39168f9 | 97 | fd0err(0x0), |
939850df | 98 | fProngID(0x0), |
99 | fSelectionMap(0) | |
b39168f9 | 100 | { |
101 | // | |
102 | // Constructor that can used for a "MC" object | |
103 | // | |
104 | ||
105 | fOwnPrimaryVtx = new AliAODVertex(vtx1); | |
106 | ||
107 | AliAODVertex *vtx = new AliAODVertex(vtx2); | |
108 | SetOwnSecondaryVtx(vtx); | |
109 | ||
110 | } | |
111 | //-------------------------------------------------------------------------- | |
3244eeed | 112 | AliAODRecoDecayHF::AliAODRecoDecayHF(const AliAODRecoDecayHF &source) : |
113 | AliAODRecoDecay(source), | |
0a65d33f | 114 | fOwnPrimaryVtx(0x0), |
460cd990 | 115 | fEventPrimaryVtx(source.fEventPrimaryVtx), |
a9b75906 | 116 | fListOfCuts(source.fListOfCuts), |
6185d025 | 117 | fd0err(0x0), |
939850df | 118 | fProngID(0x0), |
119 | fSelectionMap(source.fSelectionMap) | |
3244eeed | 120 | { |
121 | // | |
122 | // Copy constructor | |
123 | // | |
0a65d33f | 124 | if(source.GetOwnPrimaryVtx()) fOwnPrimaryVtx = new AliAODVertex(*(source.GetOwnPrimaryVtx())); |
125 | ||
3244eeed | 126 | if(source.GetNProngs()>0) { |
127 | fd0err = new Double_t[GetNProngs()]; | |
128 | memcpy(fd0err,source.fd0err,GetNProngs()*sizeof(Double_t)); | |
6185d025 | 129 | if(source.fProngID) { |
130 | fProngID = new UShort_t[GetNProngs()]; | |
131 | memcpy(fProngID,source.fProngID,GetNProngs()*sizeof(UShort_t)); | |
132 | } | |
3244eeed | 133 | } |
134 | } | |
135 | //-------------------------------------------------------------------------- | |
136 | AliAODRecoDecayHF &AliAODRecoDecayHF::operator=(const AliAODRecoDecayHF &source) | |
137 | { | |
138 | // | |
139 | // assignment operator | |
140 | // | |
141 | if(&source == this) return *this; | |
dcb444c9 | 142 | |
143 | AliAODRecoDecay::operator=(source); | |
144 | ||
460cd990 | 145 | fEventPrimaryVtx = source.fEventPrimaryVtx; |
a9b75906 | 146 | fListOfCuts = source.fListOfCuts; |
939850df | 147 | fSelectionMap = source.fSelectionMap; |
460cd990 | 148 | |
0a65d33f | 149 | if(source.GetOwnPrimaryVtx()) fOwnPrimaryVtx = new AliAODVertex(*(source.GetOwnPrimaryVtx())); |
150 | ||
3244eeed | 151 | if(source.GetNProngs()>0) { |
d95be5e9 | 152 | if(source.fd0err) { |
153 | fd0err = new Double_t[GetNProngs()]; | |
154 | memcpy(fd0err,source.fd0err,GetNProngs()*sizeof(Double_t)); | |
155 | } | |
6185d025 | 156 | if(source.fProngID) { |
157 | fProngID = new UShort_t[GetNProngs()]; | |
158 | memcpy(fProngID,source.fProngID,GetNProngs()*sizeof(UShort_t)); | |
159 | } | |
3244eeed | 160 | } |
161 | return *this; | |
162 | } | |
163 | //-------------------------------------------------------------------------- | |
164 | AliAODRecoDecayHF::~AliAODRecoDecayHF() { | |
165 | // | |
166 | // Default Destructor | |
167 | // | |
168 | if(fOwnPrimaryVtx) delete fOwnPrimaryVtx; | |
169 | if(fd0err) delete [] fd0err; | |
6185d025 | 170 | if(fProngID) delete [] fProngID; |
3244eeed | 171 | } |
172 | //--------------------------------------------------------------------------- | |
ec653946 | 173 | AliKFParticle *AliAODRecoDecayHF::ApplyVertexingKF(Int_t *iprongs,Int_t nprongs,Int_t *pdgs,Bool_t topoCostraint, Double_t bzkG, Double_t *mass) const { |
174 | // | |
175 | // Applies the KF vertexer | |
176 | // Int_t iprongs[nprongs] = indices of the prongs to be used from the vertexer | |
177 | // Int_t pdgs[nprongs] = pdgs assigned to the prongs, needed to define the AliKFParticle | |
178 | // Bool_t topoCostraint = if kTRUE, the topological constraint is applied | |
179 | // Double_t bzkG = magnetic field | |
180 | // Double_t mass[2] = {mass, sigma} for the mass constraint (if mass[0]>0 the constraint is applied). | |
181 | // | |
182 | ||
183 | AliKFParticle::SetField(bzkG); | |
184 | AliKFParticle *vertexKF=0; | |
185 | ||
186 | AliKFVertex copyKF; | |
187 | Int_t nt=0,ntcheck=0; | |
188 | ||
189 | Double_t pos[3]={0.,0.,0.}; | |
190 | if(!fOwnPrimaryVtx) { | |
191 | printf("AliAODRecoDecayHF::ApplyVertexingKF(): cannot apply because primary vertex is not found\n"); | |
192 | return vertexKF; | |
193 | } | |
194 | fOwnPrimaryVtx->GetXYZ(pos); | |
195 | Int_t contr=fOwnPrimaryVtx->GetNContributors(); | |
196 | Double_t covmatrix[6]={0.,0.,0.,0.,0.,0.}; | |
197 | fOwnPrimaryVtx->GetCovarianceMatrix(covmatrix); | |
198 | Double_t chi2=fOwnPrimaryVtx->GetChi2(); | |
199 | AliESDVertex primaryVtx2(pos,covmatrix,chi2,contr,"Vertex"); | |
200 | ||
201 | ||
202 | if(topoCostraint){ | |
203 | copyKF=AliKFVertex(primaryVtx2); | |
204 | nt=primaryVtx2.GetNContributors(); | |
205 | ntcheck=nt; | |
206 | } | |
207 | ||
208 | vertexKF = new AliKFParticle(); | |
209 | for(Int_t i= 0;i<nprongs;i++){ | |
210 | Int_t ipr=iprongs[i]; | |
211 | AliAODTrack *aodTrack = (AliAODTrack*)GetDaughter(ipr); | |
212 | if(!aodTrack) { | |
213 | printf("AliAODRecoDecayHF::ApplyVertexingKF(): no daughters available\n"); | |
214 | delete vertexKF; vertexKF=NULL; | |
215 | return vertexKF; | |
216 | } | |
217 | AliKFParticle daughterKF(*aodTrack,pdgs[i]); | |
218 | vertexKF->AddDaughter(daughterKF); | |
219 | ||
220 | if(topoCostraint && nt>0){ | |
221 | //Int_t index=(Int_t)GetProngID(ipr); | |
222 | if(!aodTrack->GetUsedForPrimVtxFit()) continue; | |
223 | copyKF -= daughterKF; | |
224 | ntcheck--; | |
225 | } | |
226 | } | |
227 | ||
228 | if(topoCostraint){ | |
229 | if(ntcheck>0) { | |
230 | copyKF += (*vertexKF); | |
231 | vertexKF->SetProductionVertex(copyKF); | |
232 | } | |
233 | } | |
234 | ||
235 | if(mass[0]>0.){ | |
236 | vertexKF->SetMassConstraint(mass[0],mass[1]); | |
237 | } | |
238 | ||
239 | return vertexKF; | |
240 | } | |
ead561f0 | 241 | //--------------------------------------------------------------------------- |
242 | AliAODVertex* AliAODRecoDecayHF::RemoveDaughtersFromPrimaryVtx(AliAODEvent *aod) { | |
243 | // | |
244 | // This method returns a primary vertex without the daughter tracks of the | |
245 | // candidate and it recalculates the impact parameters and errors. | |
246 | // | |
247 | // The output vertex is created with "new". The user has to | |
248 | // set it to the candidate with SetOwnPrimaryVtx(), unset it at the end | |
249 | // of processing with UnsetOwnPrimaryVtx() and delete it. | |
250 | // If a NULL pointer is returned, the removal failed (too few tracks left). | |
251 | // | |
252 | // For the moment, the primary vertex is recalculated from scratch without | |
253 | // the daughter tracks. | |
254 | // | |
255 | ||
256 | AliAODVertex *vtxAOD = aod->GetPrimaryVertex(); | |
257 | if(!vtxAOD) return 0; | |
258 | TString title=vtxAOD->GetTitle(); | |
259 | if(!title.Contains("VertexerTracks")) return 0; | |
260 | ||
261 | ||
262 | ||
263 | AliVertexerTracks *vertexer = new AliVertexerTracks(aod->GetMagneticField()); | |
264 | ||
265 | Int_t ndg = GetNDaughters(); | |
266 | ||
267 | vertexer->SetITSMode(); | |
4a037efc | 268 | vertexer->SetMinClusters(3); |
ead561f0 | 269 | vertexer->SetConstraintOff(); |
270 | ||
271 | if(title.Contains("WithConstraint")) { | |
272 | Float_t diamondcovxy[3]; | |
273 | aod->GetDiamondCovXY(diamondcovxy); | |
274 | Double_t pos[3]={aod->GetDiamondX(),aod->GetDiamondY(),0.}; | |
275 | Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.}; | |
276 | AliESDVertex *diamond = new AliESDVertex(pos,cov,1.,1); | |
277 | vertexer->SetVtxStart(diamond); | |
278 | delete diamond; diamond=NULL; | |
279 | } | |
280 | ||
281 | Int_t skipped[10]; | |
282 | Int_t nTrksToSkip=0,id; | |
283 | AliAODTrack *t = 0; | |
284 | for(Int_t i=0; i<ndg; i++) { | |
285 | t = (AliAODTrack*)GetDaughter(i); | |
286 | id = (Int_t)t->GetID(); | |
287 | if(id<0) continue; | |
288 | skipped[nTrksToSkip++] = id; | |
289 | } | |
290 | vertexer->SetSkipTracks(nTrksToSkip,skipped); | |
291 | AliESDVertex *vtxESDNew = vertexer->FindPrimaryVertex(aod); | |
292 | ||
293 | delete vertexer; vertexer=NULL; | |
294 | ||
295 | if(!vtxESDNew) return 0; | |
296 | if(vtxESDNew->GetNContributors()<=0) { | |
297 | delete vtxESDNew; vtxESDNew=NULL; | |
298 | return 0; | |
299 | } | |
300 | ||
301 | // convert to AliAODVertex | |
302 | Double_t pos[3],cov[6],chi2perNDF; | |
303 | vtxESDNew->GetXYZ(pos); // position | |
304 | vtxESDNew->GetCovMatrix(cov); //covariance matrix | |
305 | chi2perNDF = vtxESDNew->GetChi2toNDF(); | |
306 | delete vtxESDNew; vtxESDNew=NULL; | |
307 | ||
308 | AliAODVertex *vtxAODNew = new AliAODVertex(pos,cov,chi2perNDF); | |
309 | ||
70841bf2 | 310 | RecalculateImpPars(vtxAODNew,aod); |
311 | ||
312 | return vtxAODNew; | |
313 | } | |
314 | //----------------------------------------------------------------------------------- | |
315 | void AliAODRecoDecayHF::RecalculateImpPars(AliAODVertex *vtxAODNew,AliAODEvent* aod) { | |
316 | // | |
ead561f0 | 317 | // now recalculate the daughters impact parameters |
70841bf2 | 318 | // |
ead561f0 | 319 | Double_t dz[2],covdz[3]; |
70841bf2 | 320 | for(Int_t i=0; i<GetNDaughters(); i++) { |
321 | AliAODTrack *t = (AliAODTrack*)GetDaughter(i); | |
fb690bb8 | 322 | AliExternalTrackParam etp; etp.CopyFromVTrack(t); |
4ec1ca0f | 323 | if(etp.PropagateToDCA(vtxAODNew,aod->GetMagneticField(),3.,dz,covdz)) { |
ead561f0 | 324 | fd0[i] = dz[0]; |
325 | fd0err[i] = TMath::Sqrt(covdz[0]); | |
326 | } | |
ead561f0 | 327 | } |
328 | ||
70841bf2 | 329 | return; |
ead561f0 | 330 | } |
15c6df42 | 331 | //----------------------------------------------------------------------------------- |
332 | void AliAODRecoDecayHF::Misalign(TString misal) { | |
333 | // | |
334 | // Method to smear the impact parameter of the duaghter tracks | |
335 | // and the sec. vtx position accordinlgy | |
336 | // Useful to study effect of misalignment. | |
337 | // The starting point are parameterisations of the impact parameter resolution | |
338 | // from MC and data | |
339 | // Errors on d0 and vtx are not recalculated (to be done) | |
340 | // | |
341 | if(misal=="null")return; | |
342 | Double_t pard0rphiMC[3]={36.7,36.,1.25};// d0(pt)=[0]+[1]/(pt^[2]); in micron, conversion to cm is done later | |
343 | Double_t pard0rphimisal[3]={0,0,0}; | |
344 | Double_t pard0zMC[3]={85.,130.,0.7};// d0(pt)=[0]+[1]/(pt^[2]); in micron, conversion to cm is done later | |
345 | Double_t pard0zmisal[3]={0,0,0}; | |
346 | if(misal=="data") { | |
347 | //use this to reproduce data d0(pt) trend for pions | |
348 | pard0rphimisal[0]=37.; | |
349 | pard0rphimisal[1]=37.5; | |
350 | pard0rphimisal[2]=1.25; | |
351 | pard0zmisal[0]=96.; | |
352 | pard0zmisal[1]=131.; | |
353 | pard0zmisal[2]=0.7; | |
354 | } | |
355 | else if(misal=="resB") { | |
356 | // temporary values: asymptotic value larger by a factor 1.2 w.r.t. MC | |
357 | pard0rphimisal[0]=44.4; | |
358 | pard0rphimisal[1]=37.5; | |
359 | pard0rphimisal[2]=1.25; | |
360 | pard0zmisal[0]=115.2; | |
361 | pard0zmisal[1]=131.; | |
362 | pard0zmisal[2]=0.7; | |
363 | } | |
364 | else if(misal=="resC") { | |
365 | // temporary values: slightly larger asymptotic value, larger values at low pt | |
366 | pard0rphimisal[0]=40.; | |
367 | pard0rphimisal[1]=40.; | |
368 | pard0rphimisal[2]=1.3; | |
369 | pard0zmisal[0]=125.; | |
370 | pard0zmisal[1]=131.; | |
371 | pard0zmisal[2]=0.85; | |
372 | } | |
373 | else printf("AliAODRecoDecayHF::Misalign(): wrong misalign type specified \n"); | |
374 | ||
375 | ||
376 | AliAODVertex *evVtx=0x0,*secVtx=0x0; | |
377 | Double_t evVtxPos[3]={-9999.,-9999.,-9999.},secVtxPos[3]={-9999.,9999.,9999.}; | |
378 | if(fOwnPrimaryVtx)fOwnPrimaryVtx->GetXYZ(evVtxPos); | |
379 | else { | |
380 | evVtx=(AliAODVertex*)(fEventPrimaryVtx.GetObject()); | |
381 | evVtx->GetXYZ(evVtxPos); | |
382 | } | |
383 | secVtx=(AliAODVertex*)GetSecondaryVtx(); | |
384 | secVtx->GetXYZ(secVtxPos); | |
385 | ||
386 | TVector3 v2v1(secVtxPos[0]-evVtxPos[0],secVtxPos[1]-evVtxPos[1],0.); | |
387 | ||
388 | Double_t sigmarphinull,sigmarphimisal,sigmarphiadd; | |
389 | Double_t sigmaznull,sigmazmisal,sigmazadd; | |
390 | Double_t deltad0rphi[10],deltad0z[10]; | |
391 | ||
392 | // loop on the two prongs | |
393 | for(Int_t i=0; i<fNProngs; i++) { | |
394 | sigmarphinull = pard0rphiMC[0]+pard0rphiMC[1]/TMath::Power(PtProng(i),pard0rphiMC[2]); | |
395 | sigmarphimisal = pard0rphimisal[0]+pard0rphimisal[1]/TMath::Power(PtProng(i),pard0rphimisal[2]); | |
396 | if(sigmarphimisal>sigmarphinull) { | |
397 | sigmarphiadd = TMath::Sqrt(sigmarphimisal*sigmarphimisal- | |
398 | sigmarphinull*sigmarphinull); | |
399 | deltad0rphi[i] = gRandom->Gaus(0.,sigmarphiadd); | |
400 | } else { | |
401 | deltad0rphi[i] = 0.; | |
402 | } | |
403 | ||
404 | sigmaznull = pard0zMC[0]+pard0zMC[1]/TMath::Power(PtProng(i),pard0zMC[2]); | |
405 | sigmazmisal = pard0zmisal[0]+pard0zmisal[1]/TMath::Power(PtProng(i),pard0zmisal[2]); | |
406 | if(sigmazmisal>sigmaznull) { | |
407 | sigmazadd = TMath::Sqrt(sigmazmisal*sigmazmisal- | |
408 | sigmaznull*sigmaznull); | |
409 | deltad0z[i] = gRandom->Gaus(0.,sigmazadd); | |
410 | } else { | |
411 | deltad0z[i] = 0.; | |
412 | } | |
413 | ||
414 | TVector3 pxy(fPx[i],fPy[i],0.); | |
415 | TVector3 pxycrossv2v1=pxy.Cross(v2v1); | |
416 | if( pxycrossv2v1.Z()*fd0[i] > 0 ) { | |
417 | secVtxPos[0]+=1.e-4*deltad0rphi[i]*(-fPy[i])/PtProng(i);// e-4: conversion to cm | |
418 | secVtxPos[1]+=1.e-4*deltad0rphi[i]*(+fPx[i])/PtProng(i); | |
419 | } else { | |
420 | secVtxPos[0]+=1.e-4*deltad0rphi[i]*(+fPy[i])/PtProng(i); | |
421 | secVtxPos[1]+=1.e-4*deltad0rphi[i]*(-fPx[i])/PtProng(i); | |
422 | } | |
423 | ||
424 | // change d0rphi | |
425 | fd0[i] += 1.e-4*deltad0rphi[i]; // e-4: conversion to cm | |
426 | // change secondary vertex z | |
427 | secVtxPos[2]+=0.5e-4*deltad0z[i]; | |
428 | } | |
429 | secVtx->SetX(secVtxPos[0]); | |
430 | secVtx->SetY(secVtxPos[1]); | |
431 | secVtx->SetZ(secVtxPos[2]); | |
432 | ||
433 | return; | |
434 | } |