39042aa40604a412ab2bd4266da7c957a6f55f0c
[u/mrichter/AliRoot.git] / HLT / global / AliHLTVertexer.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, 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 #include "AliHLTVertexer.h"
18 #include "AliTracker.h"
19 #include "TMath.h"
20 #include "AliESDtrack.h"
21 #include "AliESDv0.h"
22 #include "AliESDVertex.h"
23 #include "AliESDEvent.h"
24 #include "AliKFParticle.h"
25 #include "AliKFVertex.h"
26
27 ClassImp(AliHLTVertexer)
28
29 AliHLTVertexer::AliHLTVertexer():
30   fESD(0),
31   fTrackInfos(0),
32   fPrimaryVtx(),
33   fFillVtxConstrainedTracks(1)
34 {
35 }
36
37 AliHLTVertexer::AliHLTVertexer(const AliHLTVertexer & ):
38   fESD(0),
39   fTrackInfos(0),
40   fPrimaryVtx(),
41   fFillVtxConstrainedTracks(1)
42 {
43 }
44
45 void AliHLTVertexer::SetESD( AliESDEvent *event )
46 {
47   //* Fill fTrackInfo array
48
49   delete[] fTrackInfos;
50   fESD = event;
51
52   AliKFParticle::SetField( fESD->GetMagneticField() );
53
54   Int_t nESDTracks=event->GetNumberOfTracks(); 
55   fTrackInfos = new AliESDTrackInfo[ nESDTracks ];
56
57   for (Int_t iTr=0; iTr<nESDTracks; iTr++){ 
58   
59     AliESDTrackInfo &info = fTrackInfos[iTr];
60     info.fOK = 0;
61     info.fPrimUsedFlag = 0;
62     
63     //* track quality check
64
65     AliESDtrack *pTrack = event->GetTrack(iTr);    
66     if( !pTrack  ) continue;
67     if (pTrack->GetKinkIndex(0)>0) continue;
68     if ( !( pTrack->GetStatus()&AliESDtrack::kTPCin ) ) continue;
69     
70     //* Construct KFParticle for the track
71
72     info.fParticle = AliKFParticle( *pTrack->GetInnerParam(), 211 );    
73     info.fOK = 1;
74   }
75 }
76
77 void AliHLTVertexer::FindPrimaryVertex(  )
78 {
79   //* Find event primary vertex
80
81   int nTracks = fESD->GetNumberOfTracks();
82
83   const AliKFParticle **vSelected = new const AliKFParticle*[nTracks]; //* Selected particles for vertex fit
84   Int_t *vIndex = new int [nTracks];                    //* Indices of selected particles
85   Bool_t *vFlag = new bool [nTracks];                    //* Flags returned by the vertex finder
86
87   fPrimaryVtx.Initialize();
88   fPrimaryVtx.SetBeamConstraint(fESD->GetDiamondX(),fESD->GetDiamondY(),0,
89                                 TMath::Sqrt(fESD->GetSigma2DiamondX()),TMath::Sqrt(fESD->GetSigma2DiamondY()),5.3);
90   
91   Int_t nSelected = 0;
92   for( Int_t i = 0; i<nTracks; i++){ 
93     if(!fTrackInfos[i].fOK ) continue;
94     //if( fESD->GetTrack(i)->GetTPCNcls()<60  ) continue;
95     const AliKFParticle &p = fTrackInfos[i].fParticle;
96     Double_t chi = p.GetDeviationFromVertex( fPrimaryVtx );      
97     if( chi > 3.5 ) continue;
98     vSelected[nSelected] = &(fTrackInfos[i].fParticle);
99     vIndex[nSelected] = i;
100     nSelected++;  
101   }
102   fPrimaryVtx.ConstructPrimaryVertex( vSelected, nSelected, vFlag, 3. );
103   for( Int_t i = 0; i<nSelected; i++){ 
104     if( vFlag[i] ) fTrackInfos[vIndex[i]].fPrimUsedFlag = 1;
105   }
106   
107   if( fPrimaryVtx.GetNContributors()>3 ){
108     AliESDVertex vESD( fPrimaryVtx.Parameters(), fPrimaryVtx.CovarianceMatrix(), fPrimaryVtx.GetChi2(), fPrimaryVtx.GetNContributors() );
109     fESD->SetPrimaryVertexTracks( &vESD );
110
111     // relate the tracks to vertex
112
113     if( fFillVtxConstrainedTracks ){      
114       for( Int_t i = 0; i<nTracks; i++ ){
115         if( !fTrackInfos[i].fPrimUsedFlag ) continue;   
116         fESD->GetTrack(i)->RelateToVertex( &vESD, fESD->GetMagneticField(),100. );
117       }
118     }
119
120   } else {
121     for( Int_t i = 0; i<nTracks; i++)
122       fTrackInfos[i].fPrimUsedFlag = 0;
123   }
124
125
126   delete[] vSelected;
127   delete[] vIndex;
128   delete[] vFlag;
129 }
130
131 void AliHLTVertexer::FindV0s(  )
132 {
133   //* V0 finder
134
135   int nTracks = fESD->GetNumberOfTracks();
136   //AliKFVertex primVtx( *fESD->GetPrimaryVertexTracks() );
137   AliKFVertex &primVtx = fPrimaryVtx;
138   if( primVtx.GetNContributors()<3 ) return;
139
140   bool *constrainedV0   = new bool[nTracks];
141   for( Int_t iTr = 0; iTr<nTracks; iTr++ ){ 
142     AliESDTrackInfo &info = fTrackInfos[iTr];
143     info.fPrimDeviation = info.fParticle.GetDeviationFromVertex( primVtx );
144     constrainedV0[iTr] = 0;
145   }
146   
147
148   for( Int_t iTr = 0; iTr<nTracks; iTr++ ){ //* first daughter
149
150     AliESDTrackInfo &info = fTrackInfos[iTr];
151     if( !info.fOK ) continue;    
152     if( info.fParticle.GetQ() >0 ) continue;    
153     if( info.fPrimDeviation <2.5 ) continue;
154
155     for( Int_t jTr = 0; jTr<nTracks; jTr++ ){  //* second daughter
156
157       AliESDTrackInfo &jnfo = fTrackInfos[jTr];
158       if( !jnfo.fOK ) continue;
159       if( jnfo.fParticle.GetQ() < 0 ) continue;
160       if( jnfo.fPrimDeviation <2.5 ) continue;
161    
162       //* construct V0 mother
163
164       AliKFParticle v0( info.fParticle, jnfo.fParticle );     
165
166       //* check V0 Chi^2
167       
168       if( v0.GetChi2()<0 || v0.GetChi2() > 9.*v0.GetNDF() ) continue;
169
170       //* subtruct daughters from primary vertex 
171
172       AliKFVertex primVtxCopy = primVtx;    
173        
174       if( info.fPrimUsedFlag ){ 
175         if( primVtxCopy.GetNContributors()<=2 ) continue;
176         primVtxCopy -= info.fParticle;
177       }
178       if( jnfo.fPrimUsedFlag ){
179         if( primVtxCopy.GetNContributors()<=2 ) continue;
180         primVtxCopy -= jnfo.fParticle;
181       }
182       //* Check v0 Chi^2 deviation from primary vertex 
183
184       if( v0.GetDeviationFromVertex( primVtxCopy ) >3. ) continue;
185
186       //* Add V0 to primary vertex to improve the primary vertex resolution
187
188       primVtxCopy += v0;      
189
190       //* Set production vertex for V0
191
192       v0.SetProductionVertex( primVtxCopy );
193
194       //* Check chi^2 for a case
195
196       if( v0.GetChi2()<0 || v0.GetChi2()> 9.*v0.GetNDF() ) continue;
197       
198       //* Get V0 decay length with estimated error
199
200       Double_t length, sigmaLength;
201       if( v0.GetDecayLength( length, sigmaLength ) ) continue;
202
203       //* Reject V0 if it decays too close[sigma] to the primary vertex
204
205       if( length  <3.5*sigmaLength ) continue;
206
207       //* add ESD v0 
208       
209       AliESDv0 v0ESD( *fESD->GetTrack( iTr ), iTr, *fESD->GetTrack( jTr ), jTr );  
210       fESD->AddV0( &v0ESD );
211
212       // relate the tracks to vertex
213
214       if( fFillVtxConstrainedTracks ){
215         if( constrainedV0[iTr] || constrainedV0[jTr]
216             || info.fPrimDeviation < 4. || jnfo.fPrimDeviation <4. ) continue;
217         AliESDVertex vESD(v0.Parameters(), v0.CovarianceMatrix(), v0.GetChi2(), 2);
218         fESD->GetTrack(iTr)->RelateToVertex( &vESD, fESD->GetMagneticField(),100. );
219         fESD->GetTrack(jTr)->RelateToVertex( &vESD, fESD->GetMagneticField(),100. );
220         constrainedV0[iTr] = 1;
221         constrainedV0[jTr] = 1; 
222       }
223     }
224   }
225   delete[] constrainedV0;
226 }
227