ITS/AliITSV0Finder.cxx

   1 /**************************************************************************
   2  * Copyright(c) 2007-2009, 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 /* $Id$ */
  16
  17 //-------------------------------------------------------------------------
  18 //      Implementation of the on-the-fly v0 finder for ITS tracker MI
  19 //          Origin: Marian Ivanov, CERN, Marian.Ivanov@cern.ch
  20 //          Extraction from AliITStrackerMI: Andrea Dainese
  21 //          Current support and development:
  22 //-------------------------------------------------------------------------
  23
  24 #include <TMatrixD.h>
  25 #include <TTree.h>
  26 #include <TString.h>
  27 #include <TRandom.h>
  28 #include <TTreeStream.h>
  29
  30
  31 #include "AliLog.h"
  32 #include "AliESDVertex.h"
  33 #include "AliESDEvent.h"
  34 #include "AliESDtrack.h"
  35 #include "AliV0.h"
  36 #include "AliHelix.h"
  37 #include "AliITSRecPoint.h"
  38 #include "AliITSReconstructor.h"
  39 #include "AliITStrackerMI.h"
  40 #include "AliITSV0Finder.h"
  41
  42 ClassImp(AliITSV0Finder)
  43
  44 AliITSV0Finder::AliITSV0Finder() {
  45   //Default constructor
  46 }
  47 //------------------------------------------------------------------------
  48 void AliITSV0Finder::UpdateTPCV0(const AliESDEvent *event,
  49                                    AliITStrackerMI *tracker)
  50 {
  51   //
  52   //try to update, or reject TPC  V0s
  53   //
  54   TObjArray *trackHypothesys = tracker->GetTrackHypothesys();
  55
  56   Int_t nv0s = event->GetNumberOfV0s();
  57   Int_t nitstracks = trackHypothesys->GetEntriesFast();
  58
  59   for (Int_t i=0;i<nv0s;i++){
  60     AliESDv0 * vertex = event->GetV0(i);
  61     Int_t ip = vertex->GetIndex(0);
  62     Int_t im = vertex->GetIndex(1);
  63     //
  64     TObjArray * arrayp = (ip<nitstracks) ? (TObjArray*)trackHypothesys->At(ip):0;
  65     TObjArray * arraym = (im<nitstracks) ? (TObjArray*)trackHypothesys->At(im):0;
  66     AliITStrackMI * trackp = (arrayp!=0) ? (AliITStrackMI*)arrayp->At(0):0;
  67     AliITStrackMI * trackm = (arraym!=0) ? (AliITStrackMI*)arraym->At(0):0;
  68     //
  69     //
  70     if (trackp){
  71       if (trackp->GetNumberOfClusters()+trackp->GetNDeadZone()>5.5){
  72         if (trackp->GetConstrain()&&trackp->GetChi2MIP(0)<3) vertex->SetStatus(-100);
  73         if (!trackp->GetConstrain()&&trackp->GetChi2MIP(0)<2) vertex->SetStatus(-100);
  74       }
  75     }
  76
  77     if (trackm){
  78       if (trackm->GetNumberOfClusters()+trackm->GetNDeadZone()>5.5){
  79         if (trackm->GetConstrain()&&trackm->GetChi2MIP(0)<3) vertex->SetStatus(-100);
  80         if (!trackm->GetConstrain()&&trackm->GetChi2MIP(0)<2) vertex->SetStatus(-100);
  81       }
  82     }
  83     if (vertex->GetStatus()==-100) continue;
  84     //
  85     Double_t xrp[3]; vertex->GetXYZ(xrp[0],xrp[1],xrp[2]);  //I.B.
  86     Int_t clayer = tracker->GetNearestLayer(xrp);                    //I.B.
  87     vertex->SetNBefore(clayer);        //
  88     vertex->SetChi2Before(9*clayer);   //
  89     vertex->SetNAfter(6-clayer);       //
  90     vertex->SetChi2After(0);           //
  91     //
  92     if (clayer >1 ){ // calculate chi2 before vertex
  93       Float_t chi2p = 0, chi2m=0;
  94       //
  95       if (trackp){
  96         for (Int_t ilayer=0;ilayer<clayer;ilayer++){
  97           if (trackp->GetClIndex(ilayer)>0){
  98             chi2p+=trackp->GetDy(ilayer)*trackp->GetDy(ilayer)/(trackp->GetSigmaY(ilayer)*trackp->GetSigmaY(ilayer))+
  99               trackp->GetDz(ilayer)*trackp->GetDz(ilayer)/(trackp->GetSigmaZ(ilayer)*trackp->GetSigmaZ(ilayer));
 100           }
 101           else{
 102             chi2p+=9;
 103           }
 104         }
 105       }else{
 106         chi2p = 9*clayer;
 107       }
 108       //
 109       if (trackm){
 110         for (Int_t ilayer=0;ilayer<clayer;ilayer++){
 111           if (trackm->GetClIndex(ilayer)>0){
 112             chi2m+=trackm->GetDy(ilayer)*trackm->GetDy(ilayer)/(trackm->GetSigmaY(ilayer)*trackm->GetSigmaY(ilayer))+
 113               trackm->GetDz(ilayer)*trackm->GetDz(ilayer)/(trackm->GetSigmaZ(ilayer)*trackm->GetSigmaZ(ilayer));
 114           }
 115           else{
 116             chi2m+=9;
 117           }
 118         }
 119       }else{
 120         chi2m = 9*clayer;
 121       }
 122       vertex->SetChi2Before(TMath::Min(chi2p,chi2m));
 123       if (TMath::Min(chi2p,chi2m)/Float_t(clayer)<4) vertex->SetStatus(-10);  // track exist before vertex
 124     }
 125
 126     if (clayer < 5 ){ // calculate chi2 after vertex
 127       Float_t chi2p = 0, chi2m=0;
 128       //
 129       if (trackp&&TMath::Abs(trackp->GetTgl())<1.){
 130         for (Int_t ilayer=clayer;ilayer<6;ilayer++){
 131           if (trackp->GetClIndex(ilayer)>0){
 132             chi2p+=trackp->GetDy(ilayer)*trackp->GetDy(ilayer)/(trackp->GetSigmaY(ilayer)*trackp->GetSigmaY(ilayer))+
 133               trackp->GetDz(ilayer)*trackp->GetDz(ilayer)/(trackp->GetSigmaZ(ilayer)*trackp->GetSigmaZ(ilayer));
 134           }
 135           else{
 136             chi2p+=9;
 137           }
 138         }
 139       }else{
 140         chi2p = 0;
 141       }
 142       //
 143       if (trackm&&TMath::Abs(trackm->GetTgl())<1.){
 144         for (Int_t ilayer=clayer;ilayer<6;ilayer++){
 145           if (trackm->GetClIndex(ilayer)>0){
 146             chi2m+=trackm->GetDy(ilayer)*trackm->GetDy(ilayer)/(trackm->GetSigmaY(ilayer)*trackm->GetSigmaY(ilayer))+
 147               trackm->GetDz(ilayer)*trackm->GetDz(ilayer)/(trackm->GetSigmaZ(ilayer)*trackm->GetSigmaZ(ilayer));
 148           }
 149           else{
 150             chi2m+=9;
 151           }
 152         }
 153       }else{
 154         chi2m = 0;
 155       }
 156       vertex->SetChi2After(TMath::Max(chi2p,chi2m));
 157       if (TMath::Max(chi2m,chi2p)/Float_t(6-clayer)>9) vertex->SetStatus(-20);  // track not found in ITS
 158     }
 159   }
 160   //
 161 }
 162 //------------------------------------------------------------------------
 163 void AliITSV0Finder::FindV02(AliESDEvent *event,
 164                                AliITStrackerMI *tracker)
 165 {
 166   //
 167   // V0 finder
 168   //
 169   //  Cuts on DCA -  R dependent
 170   //          max distance DCA between 2 tracks cut
 171   //          maxDist = TMath::Min(kMaxDist,kMaxDist0+pvertex->GetRr()*kMaxDist);
 172   //
 173   const Float_t kMaxDist0      = 0.1;
 174   const Float_t kMaxDist1      = 0.1;
 175   const Float_t kMaxDist       = 1;
 176   const Float_t kMinPointAngle  = 0.85;
 177   const Float_t kMinPointAngle2 = 0.99;
 178   const Float_t kMinR           = 0.5;
 179   const Float_t kMaxR           = 220;
 180   //const Float_t kCausality0Cut   = 0.19;
 181   //const Float_t kLikelihood01Cut = 0.25;
 182   //const Float_t kPointAngleCut   = 0.9996;
 183   const Float_t kCausality0Cut   = 0.19;
 184   const Float_t kLikelihood01Cut = 0.45;
 185   const Float_t kLikelihood1Cut  = 0.5;
 186   const Float_t kCombinedCut     = 0.55;
 187
 188   TObjArray *trackHypothesys = tracker->GetTrackHypothesys();
 189   TObjArray *bestHypothesys  = tracker->GetBestHypothesys();
 190   TObjArray *originalTracks  = tracker->GetOriginal();
 191   //
 192   //
 193   TTreeSRedirector &cstream = *(tracker->GetDebugStreamer());
 194   Int_t ntracks    = event->GetNumberOfTracks();
 195   Int_t nitstracks = trackHypothesys->GetEntriesFast();
 196   originalTracks->Expand(ntracks);
 197   trackHypothesys->Expand(ntracks);
 198   bestHypothesys->Expand(ntracks);
 199   //
 200   AliHelix * helixes   = new AliHelix[ntracks+2];
 201   TObjArray trackarray(ntracks+2);     //array with tracks - with vertex constrain
 202   TObjArray trackarrayc(ntracks+2);    //array of "best    tracks" - without vertex constrain
 203   TObjArray trackarrayl(ntracks+2);    //array of "longest tracks" - without vertex constrain
 204   Bool_t * forbidden   = new Bool_t [ntracks+2];
 205   Int_t   *itsmap      = new Int_t  [ntracks+2];
 206   Float_t *dist        = new Float_t[ntracks+2];
 207   Float_t *normdist0   = new Float_t[ntracks+2];
 208   Float_t *normdist1   = new Float_t[ntracks+2];
 209   Float_t *normdist    = new Float_t[ntracks+2];
 210   Float_t *norm        = new Float_t[ntracks+2];
 211   Float_t *maxr        = new Float_t[ntracks+2];
 212   Float_t *minr        = new Float_t[ntracks+2];
 213   Float_t *minPointAngle= new Float_t[ntracks+2];
 214   //
 215   AliV0 *pvertex      = new AliV0;
 216   AliITStrackMI * dummy= new AliITStrackMI;
 217   dummy->SetLabel(0);
 218   AliITStrackMI  trackat0;    //temporary track for DCA calculation
 219   //
 220   Float_t primvertex[3]={tracker->GetX(),tracker->GetY(),tracker->GetZ()};
 221   //
 222   // make ITS -  ESD map
 223   //
 224   for (Int_t itrack=0;itrack<ntracks+2;itrack++) {
 225     itsmap[itrack]        = -1;
 226     forbidden[itrack]     = kFALSE;
 227     maxr[itrack]          = kMaxR;
 228     minr[itrack]          = kMinR;
 229     minPointAngle[itrack] = kMinPointAngle;
 230   }
 231   for (Int_t itrack=0;itrack<nitstracks;itrack++){
 232     AliITStrackMI * original =   (AliITStrackMI*)(originalTracks->At(itrack));
 233     Int_t           esdindex =   original->GetESDtrack()->GetID();
 234     itsmap[esdindex]         =   itrack;
 235   }
 236   //
 237   // create ITS tracks from ESD tracks if not done before
 238   //
 239   for (Int_t itrack=0;itrack<ntracks;itrack++){
 240     if (itsmap[itrack]>=0) continue;
 241     AliITStrackMI * tpctrack = new AliITStrackMI(*(event->GetTrack(itrack)));
 242     //tpctrack->fD[0] = tpctrack->GetD(GetX(),GetY());
 243     //tpctrack->fD[1] = tpctrack->GetZat(GetX())-GetZ();
 244     tpctrack->GetDZ(tracker->GetX(),tracker->GetY(),tracker->GetZ(),tpctrack->GetDP());   //I.B.
 245     if (tpctrack->GetD(0)<20 && tpctrack->GetD(1)<20){
 246       // tracks which can reach inner part of ITS
 247       // propagate track to outer its volume - with correction for material
 248       AliITStrackerMI::CorrectForTPCtoITSDeadZoneMaterial(tpctrack);
 249     }
 250     itsmap[itrack] = nitstracks;
 251     originalTracks->AddAt(tpctrack,nitstracks);
 252     nitstracks++;
 253   }
 254   //
 255   // fill temporary arrays
 256   //
 257   for (Int_t itrack=0;itrack<ntracks;itrack++){
 258     AliESDtrack *  esdtrack = event->GetTrack(itrack);
 259     Int_t          itsindex = itsmap[itrack];
 260     AliITStrackMI *original = (AliITStrackMI*)originalTracks->At(itsmap[itrack]);
 261     if (!original) continue;
 262     AliITStrackMI *bestConst  = 0;
 263     AliITStrackMI *bestLong   = 0;
 264     AliITStrackMI *best       = 0;
 265     //
 266     //
 267     TObjArray * array    = (TObjArray*)  trackHypothesys->At(itsindex);
 268     Int_t       hentries = (array==0) ?  0 : array->GetEntriesFast();
 269     // Get best track with vertex constrain
 270     for (Int_t ih=0;ih<hentries;ih++){
 271       AliITStrackMI * trackh = (AliITStrackMI*)array->At(ih);
 272       if (!trackh->GetConstrain()) continue;
 273       if (!bestConst) bestConst = trackh;
 274       if (trackh->GetNumberOfClusters()>5.0){
 275         bestConst  = trackh;                         // full track -  with minimal chi2
 276         break;
 277       }
 278       if (trackh->GetNumberOfClusters()+trackh->GetNDeadZone()<=bestConst->GetNumberOfClusters()+bestConst->GetNDeadZone())  continue;
 279       bestConst = trackh;
 280       break;
 281     }
 282     // Get best long track without vertex constrain and best track without vertex constrain
 283     for (Int_t ih=0;ih<hentries;ih++){
 284       AliITStrackMI * trackh = (AliITStrackMI*)array->At(ih);
 285       if (trackh->GetConstrain()) continue;
 286       if (!best)     best     = trackh;
 287       if (!bestLong) bestLong = trackh;
 288       if (trackh->GetNumberOfClusters()>5.0){
 289         bestLong  = trackh;                         // full track -  with minimal chi2
 290         break;
 291       }
 292       if (trackh->GetNumberOfClusters()+trackh->GetNDeadZone()<=bestLong->GetNumberOfClusters()+bestLong->GetNDeadZone())  continue;
 293       bestLong = trackh;
 294     }
 295     if (!best) {
 296       best     = original;
 297       bestLong = original;
 298     }
 299     //I.B. trackat0 = *bestLong;
 300     new (&trackat0) AliITStrackMI(*bestLong);
 301     Double_t xx,yy,zz,alpha;
 302     if (!bestLong->GetGlobalXYZat(bestLong->GetX(),xx,yy,zz)) continue;
 303     alpha = TMath::ATan2(yy,xx);
 304     if (!trackat0.Propagate(alpha,0)) continue;
 305     // calculate normalized distances to the vertex
 306     //
 307     Float_t ptfac  = (1.+100.*TMath::Abs(trackat0.GetC()));
 308     if ( bestLong->GetNumberOfClusters()>3 ){
 309       dist[itsindex]      = trackat0.GetY();
 310       norm[itsindex]      = ptfac*TMath::Sqrt(trackat0.GetSigmaY2());
 311       normdist0[itsindex] = TMath::Abs(trackat0.GetY()/norm[itsindex]);
 312       normdist1[itsindex] = TMath::Abs((trackat0.GetZ()-primvertex[2])/(ptfac*TMath::Sqrt(trackat0.GetSigmaZ2())));
 313       normdist[itsindex]  = TMath::Sqrt(normdist0[itsindex]*normdist0[itsindex]+normdist1[itsindex]*normdist1[itsindex]);
 314       if (!bestConst){
 315         if (bestLong->GetNumberOfClusters()+bestLong->GetNDeadZone()<6) normdist[itsindex]*=2.;
 316         if (bestLong->GetNumberOfClusters()+bestLong->GetNDeadZone()<5) normdist[itsindex]*=2.;
 317         if (bestLong->GetNumberOfClusters()+bestLong->GetNDeadZone()<4) normdist[itsindex]*=2.;
 318       }else{
 319         if (bestConst->GetNumberOfClusters()+bestConst->GetNDeadZone()<6) normdist[itsindex]*=1.5;
 320         if (bestConst->GetNumberOfClusters()+bestConst->GetNDeadZone()<5) normdist[itsindex]*=1.5;
 321       }
 322     }
 323     else{
 324       if (bestConst&&bestConst->GetNumberOfClusters()+bestConst->GetNDeadZone()>4.5){
 325         dist[itsindex] = bestConst->GetD(0);
 326         norm[itsindex] = bestConst->GetDnorm(0);
 327         normdist0[itsindex] = TMath::Abs(bestConst->GetD(0)/norm[itsindex]);
 328         normdist1[itsindex] = TMath::Abs(bestConst->GetD(0)/norm[itsindex]);
 329         normdist[itsindex]  = TMath::Sqrt(normdist0[itsindex]*normdist0[itsindex]+normdist1[itsindex]*normdist1[itsindex]);
 330       }else{
 331         dist[itsindex]      = trackat0.GetY();
 332         norm[itsindex]      = ptfac*TMath::Sqrt(trackat0.GetSigmaY2());
 333         normdist0[itsindex] = TMath::Abs(trackat0.GetY()/norm[itsindex]);
 334         normdist1[itsindex] = TMath::Abs((trackat0.GetZ()-primvertex[2])/(ptfac*TMath::Sqrt(trackat0.GetSigmaZ2())));
 335         normdist[itsindex]  = TMath::Sqrt(normdist0[itsindex]*normdist0[itsindex]+normdist1[itsindex]*normdist1[itsindex]);
 336         if (TMath::Abs(trackat0.GetTgl())>1.05){
 337           if (normdist[itsindex]<3) forbidden[itsindex]=kTRUE;
 338           if (normdist[itsindex]>3) {
 339             minr[itsindex] = TMath::Max(Float_t(40.),minr[itsindex]);
 340           }
 341         }
 342       }
 343     }
 344     //
 345     //-----------------------------------------------------------
 346     //Forbid primary track candidates -
 347     //
 348     //treetr->SetAlias("forbidden0","Tr0.fN<4&&Tr1.fN+Tr1.fNDeadZone>4.5");
 349     //treetr->SetAlias("forbidden1","ND<3&&Tr1.fN+Tr1.fNDeadZone>5.5");
 350     //treetr->SetAlias("forbidden2","ND<2&&Tr1.fClIndex[0]>0&&Tr1.fClIndex[0]>0");
 351     //treetr->SetAlias("forbidden3","ND<1&&Tr1.fClIndex[0]>0");
 352     //treetr->SetAlias("forbidden4","ND<4&&Tr1.fNormChi2[0]<2");
 353     //treetr->SetAlias("forbidden5","ND<5&&Tr1.fNormChi2[0]<1");
 354     //-----------------------------------------------------------
 355     if (bestConst){
 356       if (bestLong->GetNumberOfClusters()<4       && bestConst->GetNumberOfClusters()+bestConst->GetNDeadZone()>4.5)               forbidden[itsindex]=kTRUE;
 357       if (normdist[itsindex]<3 && bestConst->GetNumberOfClusters()+bestConst->GetNDeadZone()>5.5)               forbidden[itsindex]=kTRUE;
 358       if (normdist[itsindex]<2 && bestConst->GetClIndex(0)>0 && bestConst->GetClIndex(1)>0 ) forbidden[itsindex]=kTRUE;
 359       if (normdist[itsindex]<1 && bestConst->GetClIndex(0)>0)                              forbidden[itsindex]=kTRUE;
 360       if (normdist[itsindex]<4 && bestConst->GetNormChi2(0)<2)                             forbidden[itsindex]=kTRUE;
 361       if (normdist[itsindex]<5 && bestConst->GetNormChi2(0)<1)                             forbidden[itsindex]=kTRUE;
 362       if (bestConst->GetNormChi2(0)<2.5) {
 363         minPointAngle[itsindex]= 0.9999;
 364         maxr[itsindex]         = 10;
 365       }
 366     }
 367     //
 368     //forbid daughter kink candidates
 369     //
 370     if (esdtrack->GetKinkIndex(0)>0) forbidden[itsindex] = kTRUE;
 371     Bool_t isElectron = kTRUE;
 372     Bool_t isProton   = kTRUE;
 373     Double_t pid[5];
 374     esdtrack->GetESDpid(pid);
 375     for (Int_t i=1;i<5;i++){
 376       if (pid[0]<pid[i]) isElectron= kFALSE;
 377       if (pid[4]<pid[i]) isProton= kFALSE;
 378     }
 379     if (isElectron){
 380       forbidden[itsindex]=kFALSE;
 381       normdist[itsindex]*=-1;
 382     }
 383     if (isProton){
 384       if (normdist[itsindex]>2) forbidden[itsindex]=kFALSE;
 385       normdist[itsindex]*=-1;
 386     }
 387
 388     //
 389     // Causality cuts in TPC volume
 390     //
 391     if (esdtrack->GetTPCdensity(0,10) >0.6)  maxr[itsindex] = TMath::Min(Float_t(110),maxr[itsindex]);
 392     if (esdtrack->GetTPCdensity(10,30)>0.6)  maxr[itsindex] = TMath::Min(Float_t(120),maxr[itsindex]);
 393     if (esdtrack->GetTPCdensity(20,40)>0.6)  maxr[itsindex] = TMath::Min(Float_t(130),maxr[itsindex]);
 394     if (esdtrack->GetTPCdensity(30,50)>0.6)  maxr[itsindex] = TMath::Min(Float_t(140),maxr[itsindex]);
 395     //
 396     if (esdtrack->GetTPCdensity(0,60)<0.4&&bestLong->GetNumberOfClusters()<3) minr[itsindex]=100;
 397     //
 398     //
 399     if (kFALSE){
 400       cstream<<"Track"<<
 401         "Tr0.="<<best<<
 402         "Tr1.="<<((bestConst)? bestConst:dummy)<<
 403         "Tr2.="<<bestLong<<
 404         "Tr3.="<<&trackat0<<
 405         "Esd.="<<esdtrack<<
 406         "Dist="<<dist[itsindex]<<
 407         "ND0="<<normdist0[itsindex]<<
 408         "ND1="<<normdist1[itsindex]<<
 409         "ND="<<normdist[itsindex]<<
 410         "Pz="<<primvertex[2]<<
 411         "Forbid="<<forbidden[itsindex]<<
 412         "\n";
 413       //
 414     }
 415     trackarray.AddAt(best,itsindex);
 416     trackarrayc.AddAt(bestConst,itsindex);
 417     trackarrayl.AddAt(bestLong,itsindex);
 418     new (&helixes[itsindex]) AliHelix(*best);
 419   }
 420   //
 421   //
 422   //
 423   // first iterration of V0 finder
 424   //
 425   for (Int_t iesd0=0;iesd0<ntracks;iesd0++){
 426     Int_t itrack0 = itsmap[iesd0];
 427     if (forbidden[itrack0]) continue;
 428     AliITStrackMI * btrack0 = (AliITStrackMI*)trackarray.At(itrack0);
 429     if (!btrack0) continue;
 430     if (btrack0->GetSign()>0 && !AliITSReconstructor::GetRecoParam()->GetStoreLikeSignV0s()) continue;
 431     AliITStrackMI *trackc0 = (AliITStrackMI*)trackarrayc.At(itrack0);
 432     //
 433     for (Int_t iesd1=iesd0+1;iesd1<ntracks;iesd1++){
 434       Int_t itrack1 = itsmap[iesd1];
 435       if (forbidden[itrack1]) continue;
 436
 437       AliITStrackMI * btrack1 = (AliITStrackMI*)trackarray.At(itrack1);
 438       if (!btrack1) continue;
 439       if (btrack1->GetSign()<0 && !AliITSReconstructor::GetRecoParam()->GetStoreLikeSignV0s()) continue;
 440       Bool_t isGold = kFALSE;
 441       if (TMath::Abs(TMath::Abs(btrack0->GetLabel())-TMath::Abs(btrack1->GetLabel()))==1){
 442         isGold = kTRUE;
 443       }
 444       AliITStrackMI *trackc1 = (AliITStrackMI*)trackarrayc.At(itrack1);
 445       AliHelix &h1 = helixes[itrack0];
 446       AliHelix &h2 = helixes[itrack1];
 447       //
 448       // find linear distance
 449       Double_t rmin =0;
 450       //
 451       //
 452       //
 453       Double_t phase[2][2],radius[2];
 454       Int_t  points = h1.GetRPHIintersections(h2, phase, radius);
 455       if    (points==0)  continue;
 456       Double_t delta[2]={1000000,1000000};
 457       rmin = radius[0];
 458       h1.ParabolicDCA(h2,phase[0][0],phase[0][1],radius[0],delta[0]);
 459       if (points==2){
 460         if (radius[1]<rmin) rmin = radius[1];
 461         h1.ParabolicDCA(h2,phase[1][0],phase[1][1],radius[1],delta[1]);
 462       }
 463       rmin = TMath::Sqrt(rmin);
 464       Double_t distance = 0;
 465       Double_t radiusC  = 0;
 466       Int_t    iphase   = 0;
 467       if (points==1 || delta[0]<delta[1]){
 468         distance = TMath::Sqrt(delta[0]);
 469         radiusC  = TMath::Sqrt(radius[0]);
 470       }else{
 471         distance = TMath::Sqrt(delta[1]);
 472         radiusC  = TMath::Sqrt(radius[1]);
 473         iphase=1;
 474       }
 475       if (radiusC<TMath::Max(minr[itrack0],minr[itrack1]))    continue;
 476       if (radiusC>TMath::Min(maxr[itrack0],maxr[itrack1]))     continue;
 477       Float_t maxDist  = TMath::Min(kMaxDist,Float_t(kMaxDist0+radiusC*kMaxDist1));
 478       if (distance>maxDist) continue;
 479       Float_t pointAngle = h1.GetPointAngle(h2,phase[iphase],primvertex);
 480       if (pointAngle<TMath::Max(minPointAngle[itrack0],minPointAngle[itrack1])) continue;
 481       //
 482       //
 483       //      Double_t distance = TestV0(h1,h2,pvertex,rmin);
 484       //
 485       //       if (distance>maxDist)           continue;
 486       //       if (pvertex->GetRr()<kMinR)     continue;
 487       //       if (pvertex->GetRr()>kMaxR)     continue;
 488       AliITStrackMI * track0=btrack0;
 489       AliITStrackMI * track1=btrack1;
 490       //      if (pvertex->GetRr()<3.5){
 491       if (radiusC<3.5){
 492         //use longest tracks inside the pipe
 493         track0 = (AliITStrackMI*)trackarrayl.At(itrack0);
 494         track1 = (AliITStrackMI*)trackarrayl.At(itrack1);
 495       }
 496       //
 497       //
 498       pvertex->SetParamN(*track0);
 499       pvertex->SetParamP(*track1);
 500       pvertex->Update(primvertex);
 501       pvertex->SetClusters(track0->ClIndex(),track1->ClIndex());  // register clusters
 502
 503       if (pvertex->GetRr()<kMinR) continue;
 504       if (pvertex->GetRr()>kMaxR) continue;
 505       if (pvertex->GetV0CosineOfPointingAngle()<kMinPointAngle) continue;
 506 //Bo:      if (pvertex->GetDist2()>maxDist) continue;
 507       if (pvertex->GetDcaV0Daughters()>maxDist) continue;
 508 //Bo:        pvertex->SetLab(0,track0->GetLabel());
 509 //Bo:        pvertex->SetLab(1,track1->GetLabel());
 510       pvertex->SetIndex(0,track0->GetESDtrack()->GetID());
 511       pvertex->SetIndex(1,track1->GetESDtrack()->GetID());
 512       //
 513       AliITStrackMI * htrackc0 = trackc0 ? trackc0:dummy;
 514       AliITStrackMI * htrackc1 = trackc1 ? trackc1:dummy;
 515
 516       //
 517       //
 518       TObjArray * array0b     = (TObjArray*)bestHypothesys->At(itrack0);
 519       if (!array0b&&pvertex->GetRr()<40 && TMath::Abs(track0->GetTgl())<1.1) {
 520         tracker->SetCurrentEsdTrack(itrack0);
 521         tracker->FollowProlongationTree((AliITStrackMI*)originalTracks->At(itrack0),itrack0, kFALSE);
 522       }
 523       TObjArray * array1b    = (TObjArray*)bestHypothesys->At(itrack1);
 524       if (!array1b&&pvertex->GetRr()<40 && TMath::Abs(track1->GetTgl())<1.1) {
 525         tracker->SetCurrentEsdTrack(itrack1);
 526         tracker->FollowProlongationTree((AliITStrackMI*)originalTracks->At(itrack1),itrack1, kFALSE);
 527       }
 528       //
 529       AliITStrackMI * track0b = (AliITStrackMI*)originalTracks->At(itrack0);
 530       AliITStrackMI * track1b = (AliITStrackMI*)originalTracks->At(itrack1);
 531       AliITStrackMI * track0l = (AliITStrackMI*)originalTracks->At(itrack0);
 532       AliITStrackMI * track1l = (AliITStrackMI*)originalTracks->At(itrack1);
 533
 534       Float_t minchi2before0=16;
 535       Float_t minchi2before1=16;
 536       Float_t minchi2after0 =16;
 537       Float_t minchi2after1 =16;
 538       Double_t xrp[3]; pvertex->GetXYZ(xrp[0],xrp[1],xrp[2]);  //I.B.
 539       Int_t maxLayer = tracker->GetNearestLayer(xrp);                   //I.B.
 540
 541       if (array0b) for (Int_t i=0;i<5;i++){
 542         // best track after vertex
 543         AliITStrackMI * btrack = (AliITStrackMI*)array0b->At(i);
 544         if (!btrack) continue;
 545         if (btrack->GetNumberOfClusters()>track0l->GetNumberOfClusters()) track0l = btrack;
 546         //      if (btrack->fX<pvertex->GetRr()-2.-0.5/(0.1+pvertex->GetAnglep()[2])) {
 547         if (btrack->GetX()<pvertex->GetRr()-2.) {
 548           if ( (maxLayer>i+2|| (i==0)) && btrack->GetNumberOfClusters()==(6-i)&&i<3){
 549             Float_t sumchi2= 0;
 550             Float_t sumn   = 0;
 551             if (maxLayer<3){   // take prim vertex as additional measurement
 552               if (normdist[itrack0]>0 && htrackc0){
 553                 sumchi2 += TMath::Min((3.-maxLayer)*normdist[itrack0]*normdist[itrack0],16.);
 554               }else{
 555                 sumchi2 +=  TMath::Min((3.-maxLayer)*(3*normdist[itrack0]*normdist[itrack0]+3.),16.);
 556               }
 557               sumn    +=  3-maxLayer;
 558             }
 559             for (Int_t ilayer=i;ilayer<maxLayer;ilayer++){
 560               sumn+=1.;
 561               if (!btrack->GetClIndex(ilayer)){
 562                 sumchi2+=25;
 563                 continue;
 564               }else{
 565                 Int_t c=( btrack->GetClIndex(ilayer) & 0x0fffffff);
 566                 for (Int_t itrack=0;itrack<4;itrack++){
 567                   AliITStrackerMI::AliITSlayer &layer=tracker->GetLayer(ilayer);
 568                   if (layer.GetClusterTracks(itrack,c)>=0 && layer.GetClusterTracks(itrack,c)!=itrack0){
 569                     sumchi2+=18.;  //shared cluster
 570                     break;
 571                   }
 572                 }
 573                 sumchi2+=btrack->GetDy(ilayer)*btrack->GetDy(ilayer)/(btrack->GetSigmaY(ilayer)*btrack->GetSigmaY(ilayer));
 574                 sumchi2+=btrack->GetDz(ilayer)*btrack->GetDz(ilayer)/(btrack->GetSigmaZ(ilayer)*btrack->GetSigmaZ(ilayer));
 575               }
 576             }
 577             sumchi2/=sumn;
 578             if (sumchi2<minchi2before0) minchi2before0=sumchi2;
 579           }
 580           continue;   //safety space - Geo manager will give exact layer
 581         }
 582         track0b       = btrack;
 583         minchi2after0 = btrack->GetNormChi2(i);
 584         break;
 585       }
 586       if (array1b) for (Int_t i=0;i<5;i++){
 587         // best track after vertex
 588         AliITStrackMI * btrack = (AliITStrackMI*)array1b->At(i);
 589         if (!btrack) continue;
 590         if (btrack->GetNumberOfClusters()>track1l->GetNumberOfClusters()) track1l = btrack;
 591         //      if (btrack->fX<pvertex->GetRr()-2-0.5/(0.1+pvertex->GetAnglep()[2])){
 592         if (btrack->GetX()<pvertex->GetRr()-2){
 593           if ((maxLayer>i+2 || (i==0))&&btrack->GetNumberOfClusters()==(6-i)&&(i<3)){
 594             Float_t sumchi2= 0;
 595             Float_t sumn   = 0;
 596             if (maxLayer<3){   // take prim vertex as additional measurement
 597               if (normdist[itrack1]>0 && htrackc1){
 598                 sumchi2 +=  TMath::Min((3.-maxLayer)*normdist[itrack1]*normdist[itrack1],16.);
 599               }else{
 600                 sumchi2 += TMath::Min((3.-maxLayer)*(3*normdist[itrack1]*normdist[itrack1]+3.),16.);
 601               }
 602               sumn    +=  3-maxLayer;
 603             }
 604             for (Int_t ilayer=i;ilayer<maxLayer;ilayer++){
 605               sumn+=1.;
 606               if (!btrack->GetClIndex(ilayer)){
 607                 sumchi2+=30;
 608                 continue;
 609               }else{
 610                 Int_t c=( btrack->GetClIndex(ilayer) & 0x0fffffff);
 611                 for (Int_t itrack=0;itrack<4;itrack++){
 612                   AliITStrackerMI::AliITSlayer &layer=tracker->GetLayer(ilayer);
 613                   if (layer.GetClusterTracks(itrack,c)>=0 && layer.GetClusterTracks(itrack,c)!=itrack1){
 614                     sumchi2+=18.;  //shared cluster
 615                     break;
 616                   }
 617                 }
 618                 sumchi2+=btrack->GetDy(ilayer)*btrack->GetDy(ilayer)/(btrack->GetSigmaY(ilayer)*btrack->GetSigmaY(ilayer));
 619                 sumchi2+=btrack->GetDz(ilayer)*btrack->GetDz(ilayer)/(btrack->GetSigmaZ(ilayer)*btrack->GetSigmaZ(ilayer));
 620               }
 621             }
 622             sumchi2/=sumn;
 623             if (sumchi2<minchi2before1) minchi2before1=sumchi2;
 624           }
 625           continue;   //safety space - Geo manager will give exact layer
 626         }
 627         track1b = btrack;
 628         minchi2after1 = btrack->GetNormChi2(i);
 629         break;
 630       }
 631       //
 632       // position resolution - used for DCA cut
 633       Float_t sigmad = track0b->GetSigmaY2()+track0b->GetSigmaZ2()+track1b->GetSigmaY2()+track1b->GetSigmaZ2()+
 634         (track0b->GetX()-pvertex->GetRr())*(track0b->GetX()-pvertex->GetRr())*(track0b->GetSigmaSnp2()+track0b->GetSigmaTgl2())+
 635         (track1b->GetX()-pvertex->GetRr())*(track1b->GetX()-pvertex->GetRr())*(track1b->GetSigmaSnp2()+track1b->GetSigmaTgl2());
 636       sigmad =TMath::Sqrt(sigmad)+0.04;
 637       if (pvertex->GetRr()>50){
 638         Double_t cov0[15],cov1[15];
 639         track0b->GetESDtrack()->GetInnerExternalCovariance(cov0);
 640         track1b->GetESDtrack()->GetInnerExternalCovariance(cov1);
 641         sigmad = cov0[0]+cov0[2]+cov1[0]+cov1[2]+
 642           (80.-pvertex->GetRr())*(80.-pvertex->GetRr())*(cov0[5]+cov0[9])+
 643           (80.-pvertex->GetRr())*(80.-pvertex->GetRr())*(cov1[5]+cov1[9]);
 644         sigmad =TMath::Sqrt(sigmad)+0.05;
 645       }
 646       //
 647       AliV0 vertex2;
 648       vertex2.SetParamN(*track0b);
 649       vertex2.SetParamP(*track1b);
 650       vertex2.Update(primvertex);
 651       //Bo:      if (vertex2.GetDist2()<=pvertex->GetDist2()&&(vertex2.GetV0CosineOfPointingAngle()>=pvertex->GetV0CosineOfPointingAngle())){
 652       if (vertex2.GetDcaV0Daughters()<=pvertex->GetDcaV0Daughters()&&(vertex2.GetV0CosineOfPointingAngle()>=pvertex->GetV0CosineOfPointingAngle())){
 653         pvertex->SetParamN(*track0b);
 654         pvertex->SetParamP(*track1b);
 655         pvertex->Update(primvertex);
 656         pvertex->SetClusters(track0b->ClIndex(),track1b->ClIndex());  // register clusters
 657         pvertex->SetIndex(0,track0->GetESDtrack()->GetID());
 658         pvertex->SetIndex(1,track1->GetESDtrack()->GetID());
 659       }
 660       pvertex->SetDistSigma(sigmad);
 661       //Bo:      pvertex->SetDistNorm(pvertex->GetDist2()/sigmad);
 662       pvertex->SetNormDCAPrim(normdist[itrack0],normdist[itrack1]);
 663       //
 664       // define likelihhod and causalities
 665       //
 666       Float_t pa0=1, pa1=1, pb0=0.26, pb1=0.26;
 667       if (maxLayer<1){
 668         Float_t fnorm0 = normdist[itrack0];
 669         if (fnorm0<0) fnorm0*=-3;
 670         Float_t fnorm1 = normdist[itrack1];
 671         if (fnorm1<0) fnorm1*=-3;
 672         if ((pvertex->GetAnglep()[2]>0.1) || ( (pvertex->GetRr()<10.5)&& pvertex->GetAnglep()[2]>0.05 ) || (pvertex->GetRr()<3)){
 673           pb0    =  TMath::Exp(-TMath::Min(fnorm0,Float_t(16.))/12.);
 674           pb1    =  TMath::Exp(-TMath::Min(fnorm1,Float_t(16.))/12.);
 675         }
 676         pvertex->SetChi2Before(normdist[itrack0]);
 677         pvertex->SetChi2After(normdist[itrack1]);
 678         pvertex->SetNAfter(0);
 679         pvertex->SetNBefore(0);
 680       }else{
 681         pvertex->SetChi2Before(minchi2before0);
 682         pvertex->SetChi2After(minchi2before1);
 683          if (pvertex->GetAnglep()[2]>0.1 || ( pvertex->GetRr()<10.5 && pvertex->GetAnglep()[2]>0.05) || pvertex->GetRr()<3){
 684            pb0    =  TMath::Exp(-TMath::Min(minchi2before0,Float_t(16))/12.);
 685            pb1    =  TMath::Exp(-TMath::Min(minchi2before1,Float_t(16))/12.);
 686          }
 687          pvertex->SetNAfter(maxLayer);
 688          pvertex->SetNBefore(maxLayer);
 689       }
 690       if (pvertex->GetRr()<90){
 691         pa0  *= TMath::Min(track0->GetESDtrack()->GetTPCdensity(0,60),Double_t(1.));
 692         pa1  *= TMath::Min(track1->GetESDtrack()->GetTPCdensity(0,60),Double_t(1.));
 693       }
 694       if (pvertex->GetRr()<20){
 695         pa0  *= (0.2+TMath::Exp(-TMath::Min(minchi2after0,Float_t(16))/8.))/1.2;
 696         pa1  *= (0.2+TMath::Exp(-TMath::Min(minchi2after1,Float_t(16))/8.))/1.2;
 697       }
 698       //
 699       pvertex->SetCausality(pb0,pb1,pa0,pa1);
 700       //
 701       //  Likelihood calculations  - apply cuts
 702       //
 703       Bool_t v0OK = kTRUE;
 704       Float_t p12 = pvertex->GetParamP()->GetParameter()[4]*pvertex->GetParamP()->GetParameter()[4];
 705       p12        += pvertex->GetParamN()->GetParameter()[4]*pvertex->GetParamN()->GetParameter()[4];
 706       p12         = TMath::Sqrt(p12);                             // "mean" momenta
 707       Float_t    sigmap0   = 0.0001+0.001/(0.1+pvertex->GetRr());
 708       Float_t    sigmap    = 0.5*sigmap0*(0.6+0.4*p12);           // "resolution: of point angle - as a function of radius and momenta
 709
 710       Float_t causalityA  = (1.0-pvertex->GetCausalityP()[0])*(1.0-pvertex->GetCausalityP()[1]);
 711       Float_t causalityB  = TMath::Sqrt(TMath::Min(pvertex->GetCausalityP()[2],Double_t(0.7))*
 712                                         TMath::Min(pvertex->GetCausalityP()[3],Double_t(0.7)));
 713       //
 714       //Bo:      Float_t likelihood0 = (TMath::Exp(-pvertex->GetDistNorm())+0.1) *(pvertex->GetDist2()<0.5)*(pvertex->GetDistNorm()<5);
 715       Float_t lDistNorm = pvertex->GetDcaV0Daughters()/pvertex->GetDistSigma();
 716       Float_t likelihood0 = (TMath::Exp(-lDistNorm)+0.1) *(pvertex->GetDcaV0Daughters()<0.5)*(lDistNorm<5);
 717
 718       Float_t likelihood1 = TMath::Exp(-(1.0001-pvertex->GetV0CosineOfPointingAngle())/sigmap)+
 719         0.4*TMath::Exp(-(1.0001-pvertex->GetV0CosineOfPointingAngle())/(4.*sigmap))+
 720         0.4*TMath::Exp(-(1.0001-pvertex->GetV0CosineOfPointingAngle())/(8.*sigmap))+
 721         0.1*TMath::Exp(-(1.0001-pvertex->GetV0CosineOfPointingAngle())/0.01);
 722       //
 723       if (causalityA<kCausality0Cut)                                          v0OK = kFALSE;
 724       if (TMath::Sqrt(likelihood0*likelihood1)<kLikelihood01Cut)              v0OK = kFALSE;
 725       if (likelihood1<kLikelihood1Cut)                                        v0OK = kFALSE;
 726       if (TMath::Power(likelihood0*likelihood1*causalityB,0.33)<kCombinedCut) v0OK = kFALSE;
 727
 728       //
 729       //
 730       if (kFALSE){
 731         Bool_t gold = TMath::Abs(TMath::Abs(track0->GetLabel())-TMath::Abs(track1->GetLabel()))==1;
 732         cstream<<"It0"<<
 733           "Tr0.="<<track0<<                       //best without constrain
 734           "Tr1.="<<track1<<                       //best without constrain
 735           "Tr0B.="<<track0b<<                     //best without constrain  after vertex
 736           "Tr1B.="<<track1b<<                     //best without constrain  after vertex
 737           "Tr0C.="<<htrackc0<<                    //best with constrain     if exist
 738           "Tr1C.="<<htrackc1<<                    //best with constrain     if exist
 739           "Tr0L.="<<track0l<<                     //longest best
 740           "Tr1L.="<<track1l<<                     //longest best
 741           "Esd0.="<<track0->GetESDtrack()<<           // esd track0 params
 742           "Esd1.="<<track1->GetESDtrack()<<           // esd track1 params
 743           "V0.="<<pvertex<<                       //vertex properties
 744           "V0b.="<<&vertex2<<                       //vertex properties at "best" track
 745           "ND0="<<normdist[itrack0]<<             //normalize distance for track0
 746           "ND1="<<normdist[itrack1]<<             //normalize distance for track1
 747           "Gold="<<gold<<                         //
 748           //      "RejectBase="<<rejectBase<<             //rejection in First itteration
 749           "OK="<<v0OK<<
 750           "rmin="<<rmin<<
 751           "sigmad="<<sigmad<<
 752           "\n";
 753       }
 754       //if (rejectBase) continue;
 755       //
 756       pvertex->SetStatus(0);
 757       //      if (rejectBase) {
 758       //        pvertex->SetStatus(-100);
 759       //}
 760       if (pvertex->GetV0CosineOfPointingAngle()>kMinPointAngle2) {
 761         //Bo:     pvertex->SetESDindexes(track0->GetESDtrack()->GetID(),track1->GetESDtrack()->GetID());
 762         pvertex->SetIndex(0,track0->GetESDtrack()->GetID());//Bo: consistency 0 for neg
 763         pvertex->SetIndex(1,track1->GetESDtrack()->GetID());//Bo: consistency 1 for pos
 764         if (v0OK){
 765           //      AliV0vertex vertexjuri(*track0,*track1);
 766           //      vertexjuri.SetESDindexes(track0->fESDtrack->GetID(),track1->fESDtrack->GetID());
 767           //      event->AddV0(&vertexjuri);
 768           pvertex->SetStatus(100);
 769         }
 770         pvertex->SetOnFlyStatus(kTRUE);
 771         pvertex->ChangeMassHypothesis(kK0Short);
 772         event->AddV0(pvertex);
 773       }
 774     }
 775   }
 776   //
 777   //
 778   // delete temporary arrays
 779   //
 780   delete[] forbidden;
 781   delete[] minPointAngle;
 782   delete[] maxr;
 783   delete[] minr;
 784   delete[] norm;
 785   delete[] normdist;
 786   delete[] normdist1;
 787   delete[] normdist0;
 788   delete[] dist;
 789   delete[] itsmap;
 790   delete[] helixes;
 791   delete   pvertex;
 792 }
 793 //------------------------------------------------------------------------
 794 void AliITSV0Finder::RefitV02(const AliESDEvent *event,
 795                                 AliITStrackerMI *tracker)
 796 {
 797   //
 798   //try to refit  V0s in the third path of the reconstruction
 799   //
 800   TTreeSRedirector &cstream = *(tracker->GetDebugStreamer());
 801   //
 802   Int_t  nv0s = event->GetNumberOfV0s();
 803   Float_t primvertex[3]={tracker->GetX(),tracker->GetY(),tracker->GetZ()};
 804   AliV0 v0temp;
 805   for (Int_t iv0 = 0; iv0<nv0s;iv0++){
 806     AliV0 * v0mi = (AliV0*)event->GetV0(iv0);
 807     if (!v0mi) continue;
 808     Int_t     itrack0   = v0mi->GetIndex(0);
 809     Int_t     itrack1   = v0mi->GetIndex(1);
 810     AliESDtrack *esd0   = event->GetTrack(itrack0);
 811     AliESDtrack *esd1   = event->GetTrack(itrack1);
 812     if (!esd0||!esd1) continue;
 813     AliITStrackMI tpc0(*esd0);
 814     AliITStrackMI tpc1(*esd1);
 815     Double_t x,y,z; v0mi->GetXYZ(x,y,z); //I.B.
 816     Double_t alpha =TMath::ATan2(y,x);   //I.B.
 817     if (v0mi->GetRr()>85){
 818       if (tpc0.Propagate(alpha,v0mi->GetRr())&&tpc1.Propagate(alpha,v0mi->GetRr())){
 819         v0temp.SetParamN(tpc0);
 820         v0temp.SetParamP(tpc1);
 821         v0temp.Update(primvertex);
 822         if (kFALSE) cstream<<"Refit"<<
 823           "V0.="<<v0mi<<
 824           "V0refit.="<<&v0temp<<
 825           "Tr0.="<<&tpc0<<
 826           "Tr1.="<<&tpc1<<
 827           "\n";
 828         //Bo:   if (v0temp.GetDist2()<v0mi->GetDist2() || v0temp.GetV0CosineOfPointingAngle()>v0mi->GetV0CosineOfPointingAngle()){
 829         if (v0temp.GetDcaV0Daughters()<v0mi->GetDcaV0Daughters() || v0temp.GetV0CosineOfPointingAngle()>v0mi->GetV0CosineOfPointingAngle()){
 830           v0mi->SetParamN(tpc0);
 831           v0mi->SetParamP(tpc1);
 832           v0mi->Update(primvertex);
 833         }
 834       }
 835       continue;
 836     }
 837     if (v0mi->GetRr()>35){
 838       AliITStrackerMI::CorrectForTPCtoITSDeadZoneMaterial(&tpc0);
 839       AliITStrackerMI::CorrectForTPCtoITSDeadZoneMaterial(&tpc1);
 840       if (tpc0.Propagate(alpha,v0mi->GetRr())&&tpc1.Propagate(alpha,v0mi->GetRr())){
 841         v0temp.SetParamN(tpc0);
 842         v0temp.SetParamP(tpc1);
 843         v0temp.Update(primvertex);
 844         if (kFALSE) cstream<<"Refit"<<
 845           "V0.="<<v0mi<<
 846           "V0refit.="<<&v0temp<<
 847           "Tr0.="<<&tpc0<<
 848           "Tr1.="<<&tpc1<<
 849           "\n";
 850         //Bo:   if (v0temp.GetDist2()<v0mi->GetDist2() || v0temp.GetV0CosineOfPointingAngle()>v0mi->GetV0CosineOfPointingAngle()){
 851         if (v0temp.GetDcaV0Daughters()<v0mi->GetDcaV0Daughters() || v0temp.GetV0CosineOfPointingAngle()>v0mi->GetV0CosineOfPointingAngle()){
 852           v0mi->SetParamN(tpc0);
 853           v0mi->SetParamP(tpc1);
 854           v0mi->Update(primvertex);
 855         }
 856       }
 857       continue;
 858     }
 859     AliITStrackerMI::CorrectForTPCtoITSDeadZoneMaterial(&tpc0);
 860     AliITStrackerMI::CorrectForTPCtoITSDeadZoneMaterial(&tpc1);
 861     //    if (tpc0.Propagate(alpha,v0mi->GetRr())&&tpc1.Propagate(alpha,v0mi->GetRr())){
 862     if (tracker->RefitAt(v0mi->GetRr(),&tpc0, v0mi->GetClusters(0)) &&
 863         tracker->RefitAt(v0mi->GetRr(),&tpc1, v0mi->GetClusters(1))){
 864       v0temp.SetParamN(tpc0);
 865       v0temp.SetParamP(tpc1);
 866       v0temp.Update(primvertex);
 867       if (kFALSE) cstream<<"Refit"<<
 868         "V0.="<<v0mi<<
 869         "V0refit.="<<&v0temp<<
 870         "Tr0.="<<&tpc0<<
 871         "Tr1.="<<&tpc1<<
 872         "\n";
 873       //Bo:      if (v0temp.GetDist2()<v0mi->GetDist2() || v0temp.GetV0CosineOfPointingAngle()>v0mi->GetV0CosineOfPointingAngle()){
 874       if (v0temp.GetDcaV0Daughters()<v0mi->GetDcaV0Daughters() || v0temp.GetV0CosineOfPointingAngle()>v0mi->GetV0CosineOfPointingAngle()){
 875         v0mi->SetParamN(tpc0);
 876         v0mi->SetParamP(tpc1);
 877         v0mi->Update(primvertex);
 878       }
 879     }
 880   }
 881 }
 882 //------------------------------------------------------------------------