]> git.uio.no Git - u/mrichter/AliRoot.git/blob - HLT/global/AliHLTGlobalVertexerComponent.cxx
Add method as parameter output
[u/mrichter/AliRoot.git] / HLT / global / AliHLTGlobalVertexerComponent.cxx
1 //**************************************************************************
2 //* This file is property of and copyright by the ALICE HLT Project        * 
3 //* ALICE Experiment at CERN, All rights reserved.                         *
4 //*                                                                        *
5 //* Primary Authors: S.Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de>    *
6 //*                  for The ALICE HLT Project.                            *
7 //*                                                                        *
8 //* Permission to use, copy, modify and distribute this software and its   *
9 //* documentation strictly for non-commercial purposes is hereby granted   *
10 //* without fee, provided that the above copyright notice appears in all   *
11 //* copies and that both the copyright notice and this permission notice   *
12 //* appear in the supporting documentation. The authors make no claims     *
13 //* about the suitability of this software for any purpose. It is          *
14 //* provided "as is" without express or implied warranty.                  *
15 //**************************************************************************
16
17 /** @file   AliHLTGlobalVertexerComponent.cxx
18     @author Sergey Gorbunov
19     @brief  Component for reconstruct primary vertex and V0's
20 */
21
22 #include "AliHLTGlobalVertexerComponent.h"
23 #include "AliHLTDataTypes.h"
24 #include "AliHLTExternalTrackParam.h"
25 #include "AliHLTGlobalBarrelTrack.h"
26 #include "AliCDBEntry.h"
27 #include "AliCDBManager.h"
28 #include <TFile.h>
29 #include <TString.h>
30 #include "TObjString.h"
31 #include "TObjArray.h"
32 #include "AliESDEvent.h"
33 #include "AliESDtrack.h"
34 #include "AliESDVertex.h"
35 #include "AliESDv0.h"
36 #include "AliHLTMessage.h"
37 #include "TMath.h"
38 #include "AliKFParticle.h"
39 #include "AliKFVertex.h"
40 #include "TStopwatch.h"
41
42 using namespace std;
43
44 /** ROOT macro for the implementation of ROOT specific class methods */
45 ClassImp(AliHLTGlobalVertexerComponent)
46
47 AliHLTGlobalVertexerComponent::AliHLTGlobalVertexerComponent()
48 :
49   fNTracks(0),
50   fTrackInfos(0),
51   fPrimaryVtx(),   
52   fFitTracksToVertex(1),
53   fConstrainedTrackDeviation(4.),
54   fV0DaughterPrimDeviation( 2.5 ),
55   fV0PrimDeviation( 3.5 ),
56   fV0Chi(3.5),
57   fV0DecayLengthInSigmas(3.),
58   fV0TimeLimit(10.e-3), 
59   fBenchmark("GlobalVertexer")
60 {
61   // see header file for class documentation
62   // or
63   // refer to README to build package
64   // or
65   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
66
67 }
68
69 AliHLTGlobalVertexerComponent::~AliHLTGlobalVertexerComponent()
70 {
71   // see header file for class documentation
72
73   if( fTrackInfos ) delete[] fTrackInfos;
74 }
75
76 // Public functions to implement AliHLTComponent's interface.
77 // These functions are required for the registration process
78
79 const char* AliHLTGlobalVertexerComponent::GetComponentID()
80 {
81   // see header file for class documentation
82   
83   return "GlobalVertexer";
84 }
85
86 void AliHLTGlobalVertexerComponent::GetInputDataTypes(AliHLTComponentDataTypeList& list)
87 {
88   // see header file for class documentation
89   list.clear();
90   list.push_back( kAliHLTDataTypeESDObject );
91   list.push_back( kAliHLTDataTypeTrack|kAliHLTDataOriginITS );
92   list.push_back( kAliHLTDataTypeTrack|kAliHLTDataOriginTPC );
93 }
94
95 AliHLTComponentDataType AliHLTGlobalVertexerComponent::GetOutputDataType()
96 {
97   // see header file for class documentation
98   return kAliHLTMultipleDataType;
99 }
100
101 int AliHLTGlobalVertexerComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList)
102
103 {
104   // see header file for class documentation
105   tgtList.clear();
106   tgtList.push_back( kAliHLTDataTypeGlobalVertexer|kAliHLTDataOriginOut);
107   tgtList.push_back( kAliHLTDataTypeESDVertex|kAliHLTDataOriginOut);
108   tgtList.push_back( kAliHLTDataTypeESDObject|kAliHLTDataOriginOut);
109   return tgtList.size();
110 }
111
112 void AliHLTGlobalVertexerComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
113 {
114   // see header file for class documentation
115   // XXX TODO: Find more realistic values.
116   constBase = 80000;
117   inputMultiplier = 2.;
118 }
119
120 AliHLTComponent* AliHLTGlobalVertexerComponent::Spawn()
121 {
122   // see header file for class documentation
123   return new AliHLTGlobalVertexerComponent;
124 }
125
126 int AliHLTGlobalVertexerComponent::DoInit( int argc, const char** argv )
127 {
128   // init
129
130   fBenchmark.Reset();
131   fBenchmark.SetTimer(0,"total");
132   fBenchmark.SetTimer(1,"convert");
133   fBenchmark.SetTimer(2,"vprim");
134   fBenchmark.SetTimer(3,"v0");
135   fV0TimeLimit = 10.e-3;
136
137   int iResult=0;
138   TString configuration="";
139   TString argument="";
140   for (int i=0; i<argc && iResult>=0; i++) {
141     argument=argv[i];
142     if (!configuration.IsNull()) configuration+=" ";
143     configuration+=argument;
144   }
145   
146   if (!configuration.IsNull()) {
147     iResult=Configure(configuration.Data());
148   }  
149
150   return iResult; 
151 }
152   
153 int AliHLTGlobalVertexerComponent::DoDeinit()
154 {
155   // see header file for class documentation
156
157   if( fTrackInfos ) delete[] fTrackInfos;
158
159   fTrackInfos = 0;
160   fFitTracksToVertex = 1;
161   fConstrainedTrackDeviation = 4.;
162   fV0DaughterPrimDeviation = 2.5 ;
163   fV0PrimDeviation =3.5;
164   fV0Chi = 3.5;
165   fV0DecayLengthInSigmas = 3.;
166   fV0TimeLimit = 10.e-3;
167
168   return 0;
169 }
170
171 int AliHLTGlobalVertexerComponent::DoEvent( const AliHLTComponentEventData& /*evtData*/, 
172                                             AliHLTComponentTriggerData& /*trigData*/ )
173 {
174   //cout<<"AliHLTGlobalVertexerComponent::DoEvent called"<<endl;
175  
176   if ( GetFirstInputBlock( kAliHLTDataTypeSOR ) || GetFirstInputBlock( kAliHLTDataTypeEOR ) )
177     return 0;
178
179   int iResult = 0;
180   //cout<<"SG: GlobalVertexer:: DoEvent called"<<endl;
181   fBenchmark.StartNewEvent();
182   fBenchmark.Start(0);
183   
184   vector< AliExternalTrackParam > tracks;
185   vector< int > trackId;
186   vector< pair<int,int> > v0s;
187
188   AliESDEvent *event = 0; 
189
190   for( const AliHLTComponentBlockData *i= GetFirstInputBlock(kAliHLTDataTypeESDObject); i!=NULL; i=GetNextInputBlock() ){
191     fBenchmark.AddInput(i->fSize);
192   }
193
194
195   for ( const TObject *iter = GetFirstInputObject(kAliHLTDataTypeESDObject); iter != NULL; iter = GetNextInputObject() ) {
196     event = dynamic_cast<AliESDEvent*>(const_cast<TObject*>( iter ) );
197     if( !event ) continue;
198     
199     event->GetStdContent();
200     Int_t nESDTracks=event->GetNumberOfTracks(); 
201
202     for (Int_t iTr=0; iTr<nESDTracks; iTr++){   
203       AliESDtrack *pTrack = event->GetTrack(iTr);    
204       if( !pTrack  ) continue;
205       if( pTrack->GetKinkIndex(0)>0) continue;
206       if( !( pTrack->GetStatus()&AliESDtrack::kTPCin ) ) continue;
207       tracks.push_back(*pTrack);
208       trackId.push_back(iTr);
209     }      
210     break;
211   }
212
213   if( tracks.size()==0 ){
214     
215     for (const AliHLTComponentBlockData* pBlock=GetFirstInputBlock(kAliHLTDataTypeTrack|kAliHLTDataOriginITS);
216          pBlock!=NULL; pBlock=GetNextInputBlock()) {
217
218       fBenchmark.AddInput(pBlock->fSize);
219
220       AliHLTTracksData* dataPtr = reinterpret_cast<AliHLTTracksData*>( pBlock->fPtr );
221       int nTracks = dataPtr->fCount;
222
223       AliHLTExternalTrackParam* currOutTrack = dataPtr->fTracklets;
224       for( int itr=0; itr<nTracks; itr++ ){
225         AliHLTGlobalBarrelTrack t(*currOutTrack);
226         tracks.push_back( t );
227         trackId.push_back( currOutTrack->fTrackID );
228         unsigned int dSize = sizeof( AliHLTExternalTrackParam ) + currOutTrack->fNPoints * sizeof( unsigned int );
229         currOutTrack = ( AliHLTExternalTrackParam* )( (( Byte_t * )currOutTrack) + dSize );
230       }
231     }
232   }
233
234   if( tracks.size()==0 ){
235     for (const AliHLTComponentBlockData* pBlock=GetFirstInputBlock(kAliHLTDataTypeTrack|kAliHLTDataOriginTPC);
236          pBlock!=NULL; pBlock=GetNextInputBlock()) {
237
238       fBenchmark.AddInput(pBlock->fSize);
239
240       AliHLTTracksData* dataPtr = reinterpret_cast<AliHLTTracksData*>( pBlock->fPtr );
241       int nTracks = dataPtr->fCount;      
242       AliHLTExternalTrackParam* currOutTrack = dataPtr->fTracklets;
243       for( int itr=0; itr<nTracks; itr++ ){
244         AliHLTGlobalBarrelTrack t(*currOutTrack);
245         tracks.push_back( t );
246         trackId.push_back( currOutTrack->fTrackID );
247         unsigned int dSize = sizeof( AliHLTExternalTrackParam ) + currOutTrack->fNPoints * sizeof( unsigned int );
248         currOutTrack = ( AliHLTExternalTrackParam* )( (( Byte_t * )currOutTrack) + dSize );
249       }
250     }
251   }
252
253   
254   // primary vertex & V0's 
255           
256   AliKFParticle::SetField( GetBz() );
257
258   fBenchmark.Start(1);
259
260   {  //* Fill fTrackInfo array
261
262     if( fTrackInfos ) delete[] fTrackInfos;
263     fTrackInfos = 0;
264     fNTracks=tracks.size(); 
265     fTrackInfos = new AliESDTrackInfo[ fNTracks ];
266     for (Int_t iTr=0; iTr<fNTracks; iTr++){       
267       AliESDTrackInfo &info = fTrackInfos[iTr];
268       info.fOK = 1;
269       info.fPrimUsedFlag = 0;
270       info.fParticle = AliKFParticle( tracks[iTr], 211 );
271       for( int i=0; i<8; i++ ) if( !finite(info.fParticle.GetParameter(i)) ) info.fOK = 0;
272       for( int i=0; i<36; i++ ) if( !finite(info.fParticle.GetCovariance(i)) ) info.fOK = 0;
273     }
274   }
275
276   fBenchmark.Stop(1);
277   fBenchmark.Start(2);
278   FindPrimaryVertex();
279   fBenchmark.Stop(2);
280   fBenchmark.Start(3);
281   FindV0s( v0s );
282   fBenchmark.Stop(3);
283
284   int *buf = new int[sizeof(AliHLTGlobalVertexerData)/sizeof(int)+1 + fNTracks + 2*v0s.size()];
285   AliHLTGlobalVertexerData *data = reinterpret_cast<AliHLTGlobalVertexerData*>(buf);
286
287   if( data) {  // fill the output structure
288         
289     data->fFitTracksToVertex = fFitTracksToVertex;
290     for( int i=0; i<3; i++ ) data->fPrimP[i] = fPrimaryVtx.Parameters()[i];
291     for( int i=0; i<6; i++ ) data->fPrimC[i] = fPrimaryVtx.CovarianceMatrix()[i];
292     data->fPrimChi2 = fPrimaryVtx.GetChi2();
293     data->fPrimNContributors = fPrimaryVtx.GetNContributors();
294     data->fNPrimTracks = 0;
295     for( Int_t i = 0; i<fNTracks; i++ ){
296       if( !fTrackInfos[i].fPrimUsedFlag ) continue;       
297       if( fTrackInfos[i].fPrimDeviation > fConstrainedTrackDeviation ) continue;
298       data->fTrackIndices[ (data->fNPrimTracks)++ ] =  trackId[i];
299     }
300     int *listV0 = data->fTrackIndices + data->fNPrimTracks;
301     data->fNV0s = v0s.size();
302     for( int i=0; i<data->fNV0s; i++ ){
303       listV0[2*i] = trackId[v0s[i].first];
304       listV0[2*i+1] = trackId[v0s[i].second];
305     }
306
307     unsigned int blockSize = sizeof(AliHLTGlobalVertexerData) + (data->fNPrimTracks + 2*data->fNV0s)*sizeof(int);
308
309     iResult = PushBack( reinterpret_cast<void*>(data), blockSize, kAliHLTDataTypeGlobalVertexer|kAliHLTDataOriginOut );  
310     fBenchmark.AddOutput(blockSize);
311   }  
312   
313   
314   // output the vertex if found
315   {
316     if( iResult==0 && data && data->fPrimNContributors >=3 ){
317       AliESDVertex vESD( data->fPrimP, data->fPrimC, data->fPrimChi2, data->fPrimNContributors );
318       iResult = PushBack( (TObject*) &vESD, kAliHLTDataTypeESDVertex|kAliHLTDataOriginOut,0 );
319       fBenchmark.AddOutput(GetLastObjectSize());
320     }
321   }
322
323   // output the ESD event
324   if( iResult==0 && event && data ){  
325     FillESD( event, data ); 
326     iResult = PushBack( event, kAliHLTDataTypeESDObject|kAliHLTDataOriginOut, 0);
327     fBenchmark.AddOutput(GetLastObjectSize());
328   }
329
330   delete[] buf;
331
332   fBenchmark.Stop(0);
333   HLTInfo(fBenchmark.GetStatistics());
334   return iResult;
335 }
336
337
338 int AliHLTGlobalVertexerComponent::Configure(const char* arguments)
339 {
340   
341   int iResult=0;
342   if (!arguments) return iResult;
343   
344   TString allArgs=arguments;
345   TString argument;
346   
347   TObjArray* pTokens=allArgs.Tokenize(" ");
348   int bMissingParam=0;
349
350   if (pTokens) {
351     for (int i=0; i<pTokens->GetEntries() && iResult>=0; i++) {
352       argument=((TObjString*)pTokens->At(i))->GetString();
353       if (argument.IsNull()) continue;
354       
355       if (argument.CompareTo("-fitTracksToVertex")==0) {
356         if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
357         HLTInfo("fitTracksToVertex is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
358         fFitTracksToVertex=((TObjString*)pTokens->At(i))->GetString().Atoi();
359         continue;
360       }
361       else if (argument.CompareTo("-constrainedTrackDeviation")==0) {
362         if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
363         HLTInfo("constrainedTrackDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
364         fConstrainedTrackDeviation=((TObjString*)pTokens->At(i))->GetString().Atof();
365         continue;
366       }
367       else if (argument.CompareTo("-v0DaughterPrimDeviation")==0) {
368         if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
369         HLTInfo("v0DaughterPrimDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
370         fV0DaughterPrimDeviation=((TObjString*)pTokens->At(i))->GetString().Atof();
371         continue;
372       }
373       else if (argument.CompareTo("-v0PrimDeviation")==0) {
374         if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
375         HLTInfo("v0PrimDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
376         fV0PrimDeviation=((TObjString*)pTokens->At(i))->GetString().Atof();
377         continue;
378       }
379       else if (argument.CompareTo("-v0Chi")==0) {
380         if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
381         HLTInfo("v0Chi is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
382         fV0Chi=((TObjString*)pTokens->At(i))->GetString().Atof();
383         continue;
384       }
385       else if (argument.CompareTo("-v0DecayLengthInSigmas")==0) {
386         if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
387         HLTInfo("v0DecayLengthInSigmas is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
388         fV0DecayLengthInSigmas=((TObjString*)pTokens->At(i))->GetString().Atof();
389         continue;
390       }
391       else if (argument.CompareTo("-v0MinEventRate")==0) {
392         if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
393         HLTInfo("Minimum event rate for V0 finder is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data());
394         Double_t rate = ((TObjString*)pTokens->At(i))->GetString().Atof();
395         fV0TimeLimit = (rate >0 ) ?1./rate :60; // 1 minute maximum time
396         continue;
397       }
398       else {
399         HLTError("unknown argument %s", argument.Data());
400         iResult=-EINVAL;
401         break;
402       }
403     }
404     delete pTokens;
405   }
406   if (bMissingParam) {
407     HLTError("missing parameter for argument %s", argument.Data());
408     iResult=-EINVAL;
409   }  
410   
411   return iResult;
412 }
413
414 int AliHLTGlobalVertexerComponent::Reconfigure(const char* /*cdbEntry*/, const char* /*chainId*/)
415 {
416   // see header file for class documentation
417
418   return 0; // no CDB path is set so far
419   /*
420   int iResult=0;  
421   const char* path="HLT/ConfigTPC/KryptonHistoComponent";
422   const char* defaultNotify="";
423   if (cdbEntry) {
424     path=cdbEntry;
425     defaultNotify=" (default)";
426   }
427   if (path) {
428     HLTInfo("reconfigure from entry %s%s, chain id %s", path, defaultNotify,(chainId!=NULL && chainId[0]!=0)?chainId:"<none>");
429     AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(path);
430     if (pEntry) {
431       TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
432       if (pString) {
433         HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data());
434         iResult=Configure(pString->GetString().Data());
435       } else {
436         HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
437       }
438     } else {
439       HLTError("can not fetch object \"%s\" from CDB", path);
440     }
441   }
442
443   return iResult;
444 */
445 }
446
447
448 struct AliHLTGlobalVertexerDeviation
449 {
450   int fI; // track index
451   float fD; // deviation from vertex
452
453   bool operator<(const AliHLTGlobalVertexerDeviation &a) const { return fD<a.fD; }
454 };
455
456
457 void AliHLTGlobalVertexerComponent::FindPrimaryVertex()
458 {
459   //* Find event primary vertex
460
461   fPrimaryVtx.Initialize();
462   //fPrimaryVtx.SetBeamConstraint(fESD->GetDiamondX(),fESD->GetDiamondY(),0,
463   //TMath::Sqrt(fESD->GetSigma2DiamondX()),TMath::Sqrt(fESD->GetSigma2DiamondY()),5.3);
464
465   // select rough region (in sigmas) in which the vertex could be found, all tracks outside these limits are rejected
466   // from the primary vertex finding
467   fPrimaryVtx.SetBeamConstraint( 0, 0, 0, 3., 3., 5.3 );
468
469   const AliKFParticle **vSelected = new const AliKFParticle*[fNTracks]; //* Selected particles for the vertex fit
470   AliHLTGlobalVertexerDeviation *dev = new AliHLTGlobalVertexerDeviation[fNTracks];  
471
472   Int_t nSelected = 0;
473   
474   for( Int_t i = 0; i<fNTracks; i++){ 
475     if(!fTrackInfos[i].fOK ) continue;
476     //if( fESD->GetTrack(i)->GetTPCNcls()<60  ) continue;
477     const AliKFParticle &p = fTrackInfos[i].fParticle;
478     Double_t chi = p.GetDeviationFromVertex( fPrimaryVtx );      
479     if( chi > fConstrainedTrackDeviation ) continue;
480     dev[nSelected].fI = i; 
481     dev[nSelected].fD = chi; 
482     vSelected[nSelected] = &(fTrackInfos[i].fParticle);
483     nSelected++;  
484   }    
485
486   // fit
487
488   while( nSelected>2 ){ 
489
490     //* Primary vertex finder with rejection of outliers
491
492     for( Int_t i = 0; i<nSelected; i++){ 
493       vSelected[i] = &(fTrackInfos[dev[i].fI].fParticle);      
494     }
495
496     double xv = fPrimaryVtx.GetX();
497     double yv = fPrimaryVtx.GetY();
498     double zv = fPrimaryVtx.GetZ(); // values from previous iteration of calculations          
499     fPrimaryVtx.Initialize();
500     fPrimaryVtx.SetBeamConstraint( 0, 0, 0, 3., 3., 5.3 );
501     fPrimaryVtx.SetVtxGuess( xv, yv, zv );    
502
503     fPrimaryVtx.Construct( vSelected, nSelected, 0, -1, 1 ); // refilled for every iteration
504     
505     for( Int_t it=0; it<nSelected; it++ ){ 
506       const AliKFParticle &p = fTrackInfos[dev[it].fI].fParticle;
507       if( nSelected <= 20 ){
508         AliKFVertex tmp =  fPrimaryVtx - p; // exclude the current track from the sample and recalculate the vertex
509         dev[it].fD = p.GetDeviationFromVertex( tmp );
510       } else {
511         dev[it].fD = p.GetDeviationFromVertex( fPrimaryVtx );   
512       }
513     }
514     sort(dev,dev+nSelected); // sort tracks with increasing chi2 (used for rejection)
515        
516     int nRemove = (int) ( 0.3*nSelected );  //remove 30% of the tracks (done for performance, only if there are more than 20 tracks)  
517     if( nSelected - nRemove <=20 ) nRemove = 1;  // removal based on the chi2 of every track   
518     int firstRemove = nSelected - nRemove;
519     while( firstRemove<nSelected ){
520       if( dev[firstRemove].fD >= fConstrainedTrackDeviation ) break;
521       firstRemove++;
522     }
523     if( firstRemove>=nSelected ) break;
524     nSelected = firstRemove;
525   }
526
527   for( Int_t i = 0; i<fNTracks; i++){ 
528     fTrackInfos[i].fPrimUsedFlag = 0;
529   }
530
531   if( nSelected < 3 ){  // no vertex for fewer than 3 contributors  
532     fPrimaryVtx.NDF() = -3;
533     fPrimaryVtx.Chi2() = 0;
534     nSelected = 0;
535   }
536   
537   for( Int_t i = 0; i<nSelected; i++){ 
538     AliESDTrackInfo &info = fTrackInfos[dev[i].fI];
539     info.fPrimUsedFlag = 1;
540     info.fPrimDeviation = dev[i].fD;
541   }
542
543   for( Int_t i = 0; i<fNTracks; i++ ){
544     AliESDTrackInfo &info = fTrackInfos[i];
545     if( info.fPrimUsedFlag ) continue;
546     info.fPrimDeviation = info.fParticle.GetDeviationFromVertex( fPrimaryVtx );   
547   }
548
549   delete[] vSelected;
550   delete[] dev;
551 }
552
553
554
555 void AliHLTGlobalVertexerComponent::FindV0s( vector<pair<int,int> > &v0s  )
556 {
557   //* V0 finder
558
559   AliKFVertex &primVtx = fPrimaryVtx;
560   if( primVtx.GetNContributors()<3 ) return;
561
562   TStopwatch timer;
563   Int_t statN = 0;
564   Bool_t run = 1;
565
566   for( Int_t iTr = 0; iTr<fNTracks && run; iTr++ ){ //* first daughter
567
568     AliESDTrackInfo &info = fTrackInfos[iTr];
569     if( !info.fOK ) continue;    
570     if( info.fParticle.GetQ() >0 ) continue;    
571     if( info.fPrimDeviation < fV0DaughterPrimDeviation ) continue;
572
573     for( Int_t jTr = 0; jTr<fNTracks; jTr++ ){  //* second daughter
574       
575       
576       AliESDTrackInfo &jnfo = fTrackInfos[jTr];
577       if( !jnfo.fOK ) continue;
578       if( jnfo.fParticle.GetQ() < 0 ) continue;
579       if( jnfo.fPrimDeviation < fV0DaughterPrimDeviation ) continue;
580
581       // check the time once a while...
582
583       if( (++statN)%100 ==0 ){ 
584         if( timer.RealTime()>= fV0TimeLimit ){  run = 0; break; }
585         timer.Continue();
586       }
587
588       //* check if the particles fit
589
590       if( info.fParticle.GetDeviationFromParticle(jnfo.fParticle) > fV0Chi ) continue;
591
592       //* construct V0 mother
593
594       AliKFParticle v0( info.fParticle, jnfo.fParticle );     
595
596       //* check V0 Chi^2
597       
598       if( v0.GetChi2()<0 || v0.GetChi2() > fV0Chi*fV0Chi*v0.GetNDF() ) continue;
599
600       //* subtruct daughters from primary vertex 
601
602       AliKFVertex primVtxCopy = primVtx;    
603        
604       if( info.fPrimUsedFlag ){ 
605         if( primVtxCopy.GetNContributors()<=2 ) continue;
606         primVtxCopy -= info.fParticle;
607       }
608       if( jnfo.fPrimUsedFlag ){
609         if( primVtxCopy.GetNContributors()<=2 ) continue;
610         primVtxCopy -= jnfo.fParticle;
611       }
612       //* Check v0 Chi^2 deviation from primary vertex 
613
614       if( v0.GetDeviationFromVertex( primVtxCopy ) > fV0PrimDeviation ) continue;
615
616       //* Add V0 to primary vertex to improve the primary vertex resolution
617
618       primVtxCopy += v0;      
619
620       //* Set production vertex for V0
621
622       v0.SetProductionVertex( primVtxCopy );
623
624       //* Get V0 decay length with estimated error
625
626       Double_t length, sigmaLength;
627       if( v0.GetDecayLength( length, sigmaLength ) ) continue;
628
629       //* Reject V0 if it decays too close[sigma] to the primary vertex
630
631       if( length  < fV0DecayLengthInSigmas*sigmaLength ) continue;
632       
633       //* keep v0 
634
635       v0s.push_back(pair<int,int>(iTr,jTr));
636     }
637   }
638 }
639
640
641
642
643 void AliHLTGlobalVertexerComponent::FillESD( AliESDEvent *event, AliHLTGlobalVertexerData *data
644 )
645 {
646   //* put output of a vertexer to the esd event
647
648   Int_t nESDTracks = event->GetNumberOfTracks();
649
650   const int *listPrim = data->fTrackIndices;
651   const int *listV0 = data->fTrackIndices + data->fNPrimTracks;
652
653   std::map<int,int> mapId;
654   bool *constrainedToVtx   = new bool[nESDTracks];
655
656   for( int i=0; i<nESDTracks; i++ ){
657     constrainedToVtx[i] = 0;
658     if( !event->GetTrack(i) ) continue;
659     mapId[ event->GetTrack(i)->GetID() ] = i;
660   }
661
662   if( data->fPrimNContributors >=3 ){
663
664     AliESDVertex vESD( data->fPrimP, data->fPrimC, data->fPrimChi2, data->fPrimNContributors );
665     event->SetPrimaryVertexTPC( &vESD );
666     event->SetPrimaryVertexTracks( &vESD );
667
668     // relate tracks to the primary vertex
669
670     if( data->fFitTracksToVertex ){
671       for( Int_t i = 0; i<data->fNPrimTracks; i++ ){
672         Int_t id = listPrim[ i ];
673         map<int,int>::iterator it = mapId.find(id);
674         if( it==mapId.end() ) continue;
675         Int_t itr = it->second;
676         event->GetTrack(itr)->RelateToVertex( &vESD, event->GetMagneticField(),100. );
677         constrainedToVtx[ itr ] = 1;
678       }
679     }
680   }
681
682   //* add ESD v0s and relate tracks to v0s
683  
684
685   for( int i=0; i<data->fNV0s; i++ ){
686
687     Int_t id1 = listV0[ 2*i ];
688     Int_t id2 = listV0[ 2*i + 1];
689     map<int,int>::iterator it = mapId.find(id1);
690     if( it==mapId.end() ) continue;
691     Int_t iTr = it->second;
692     it = mapId.find(id2);
693     if( it==mapId.end() ) continue;
694     Int_t jTr = it->second;
695    
696     AliESDv0 v0( *event->GetTrack( iTr ), iTr, *event->GetTrack( jTr ), jTr );  
697     event->AddV0( &v0 );
698
699     // relate the tracks to the vertex  
700
701     if( data->fFitTracksToVertex ){
702       if( constrainedToVtx[iTr] || constrainedToVtx[jTr] ) continue;
703       double pos[3];
704       double sigma[3] = {.1,.1,.1};
705       v0.XvYvZv(pos);
706       AliESDVertex vESD(pos, sigma);
707       event->GetTrack(iTr)->RelateToVertex( &vESD, event->GetMagneticField(),100. );
708       event->GetTrack(jTr)->RelateToVertex( &vESD, event->GetMagneticField(),100. );
709       constrainedToVtx[iTr] = 1;
710       constrainedToVtx[jTr] = 1;    
711     }
712   }
713
714   delete[] constrainedToVtx;
715 }