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