X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=ITS%2FAliITStrackerV2.cxx;h=4f9955d85a9e430093d87b24b0263cce12cef405;hb=cfe729e0df60b6ba4d6d81fdd4169801d6d3c235;hp=f27b2f8deae01982b71b910fbb36ec8210807d6b;hpb=880f41b9b7f55ee837a38631a3dda0729a4686c8;p=u%2Fmrichter%2FAliRoot.git diff --git a/ITS/AliITStrackerV2.cxx b/ITS/AliITStrackerV2.cxx index f27b2f8deae..4f9955d85a9 100644 --- a/ITS/AliITStrackerV2.cxx +++ b/ITS/AliITStrackerV2.cxx @@ -15,11 +15,14 @@ //------------------------------------------------------------------------- // Implementation of the ITS tracker class -// +// It reads AliITSclusterV2 clusters and creates AliITStrackV2 tracks +// and fills with them the ESD // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch // dEdx analysis by: Boris Batyunya, JINR, Boris.Batiounia@cern.ch //------------------------------------------------------------------------- +#include + #include #include #include @@ -27,12 +30,13 @@ #include "AliITSgeom.h" #include "AliITSRecPoint.h" #include "AliTPCtrack.h" +#include "AliESD.h" #include "AliITSclusterV2.h" #include "AliITStrackerV2.h" ClassImp(AliITStrackerV2) -AliITStrackerV2::AliITSlayer AliITStrackerV2::fLayers[kMaxLayer]; // ITS layers +AliITStrackerV2::AliITSlayer AliITStrackerV2::fgLayers[kMaxLayer]; // ITS layers AliITStrackerV2::AliITStrackerV2(const AliITSgeom *geom) : AliTracker() { //-------------------------------------------------------------------- @@ -59,18 +63,20 @@ AliITStrackerV2::AliITStrackerV2(const AliITSgeom *geom) : AliTracker() { r += TMath::Sqrt(x*x + y*y); r*=0.25; - new (fLayers+i-1) AliITSlayer(r,poff,zoff,nlad,ndet); + new (fgLayers+i-1) AliITSlayer(r,poff,zoff,nlad,ndet); for (Int_t j=1; jGetTrans(i,j,k,x,y,zshift); Double_t rot[9]; g->GetRotMatrix(i,j,k,rot); - Double_t r =-x*rot[1] + y*rot[0]; if (i==1) r=-r; - Double_t phi=TMath::ATan2(rot[1],rot[0]); if (i==1) phi-=3.1415927; - phi+=0.5*TMath::Pi(); if (phi<0) phi += 2*TMath::Pi(); - AliITSdetector &det=fLayers[i-1].GetDetector((j-1)*ndet + k-1); + Double_t phi=TMath::ATan2(rot[1],rot[0])+TMath::Pi(); + phi+=TMath::Pi()/2; + if (i==1) phi+=TMath::Pi(); + Double_t cp=TMath::Cos(phi), sp=TMath::Sin(phi); + Double_t r=x*cp+y*sp; + AliITSdetector &det=fgLayers[i-1].GetDetector((j-1)*ndet + k-1); new(&det) AliITSdetector(r,phi); } } @@ -84,22 +90,26 @@ AliITStrackerV2::AliITStrackerV2(const AliITSgeom *geom) : AliTracker() { Double_t xyz[]={kXV,kYV,kZV}, ers[]={kSigmaXV,kSigmaYV,kSigmaZV}; SetVertex(xyz,ers); + + for (Int_t i=0; iGet(cname); - if (!cTree) { - Error("LoadClusters"," can't get cTree !\n"); - return 1; - } TBranch *branch=cTree->GetBranch("Clusters"); if (!branch) { - Error("LoadClusters"," can't get Clusters branch !\n"); + Error("LoadClusters"," can't get the branch !\n"); return 1; } @@ -107,21 +117,37 @@ Int_t AliITStrackerV2::LoadClusters() { branch->SetAddress(&clusters); Int_t j=0; + Int_t detector=0; for (Int_t i=0; iGetEvent(j)) continue; Int_t ncl=clusters->GetEntriesFast(); while (ncl--) { AliITSclusterV2 *c=(AliITSclusterV2*)clusters->UncheckedAt(ncl); - fLayers[i].InsertCluster(new AliITSclusterV2(*c)); + detector = c->GetDetectorIndex(); + fgLayers[i].InsertCluster(new AliITSclusterV2(*c)); } clusters->Delete(); + //add dead zone virtual "cluster" + + if (i<1){ + for (Float_t ydead = -2.; ydead < 2. ; ydead+=0.05){ + Int_t lab[4] = {0,0,0,detector}; + Int_t info[3] = {0,0,0}; + Float_t hit[5]={ydead,0,1,0.01,0}; + fgLayers[i].InsertCluster(new AliITSclusterV2(lab, hit, info)); + // + hit[1]=-7.; + fgLayers[i].InsertCluster(new AliITSclusterV2(lab, hit, info)); + } + } + } - fLayers[i].ResetRoad(); //road defined by the cluster density + // + fgLayers[i].ResetRoad(); //road defined by the cluster density } - delete cTree; //Thanks to Mariana Bondila return 0; } @@ -130,7 +156,7 @@ void AliITStrackerV2::UnloadClusters() { //-------------------------------------------------------------------- //This function unloads ITS clusters //-------------------------------------------------------------------- - for (Int_t i=0; iGetGlobalXYZat(rr,x,y,z); //if (TMath::Abs(y)PropagateTo(rr,-dr,x0r); if (!t->PropagateTo(rcd,-dcd,x0cd)) return 1; - if (!t->PropagateTo(riw,-diw,x0iw)) return 1; + if (!t->PropagateTo(riw+0.001,-diw,x0iw)) return 1; } else { ::Error("CorrectForDeadZoneMaterial","track is already in the dead zone !"); return 1; @@ -167,38 +193,166 @@ static Int_t CorrectForDeadZoneMaterial(AliITStrackV2 *t) { return 0; } -Int_t AliITStrackerV2::Clusters2Tracks(const TFile *inp, TFile *out) { +Int_t AliITStrackerV2::Clusters2Tracks(AliESD *event) { //-------------------------------------------------------------------- - //This functions reconstructs ITS tracks + // This functions reconstructs ITS tracks + // The clusters must be already loaded ! //-------------------------------------------------------------------- - TFile *in=(TFile*)inp; - TDirectory *savedir=gDirectory; + TObjArray itsTracks(15000); - if (LoadClusters()!=0) return 1; + {/* Read ESD tracks */ + Int_t nentr=event->GetNumberOfTracks(); + Info("Clusters2Tracks", "Number of ESD tracks: %d\n", nentr); + while (nentr--) { + AliESDtrack *esd=event->GetTrack(nentr); + + if ((esd->GetStatus()&AliESDtrack::kTPCin)==0) continue; + if (esd->GetStatus()&AliESDtrack::kTPCout) continue; + if (esd->GetStatus()&AliESDtrack::kITSin) continue; + + AliITStrackV2 *t=0; + try { + t=new AliITStrackV2(*esd); + } catch (const Char_t *msg) { + Warning("Clusters2Tracks",msg); + delete t; + continue; + } + if (TMath::Abs(t->GetD())>5) { + delete t; + continue; + } - if (!in->IsOpen()) { - Error("Clusters2Tracks","file with TPC tracks is not open !\n"); - return 1; + if (CorrectForDeadZoneMaterial(t)!=0) { + Warning("Clusters2Tracks", + "failed to correct for the material in the dead zone !\n"); + delete t; + continue; + } + t->fReconstructed = kFALSE; + itsTracks.AddLast(t); + } + } /* End Read ESD tracks */ + + itsTracks.Sort(); + Int_t nentr=itsTracks.GetEntriesFast(); + fTrackHypothesys.Expand(nentr); + Int_t ntrk=0; + for (fPass=0; fPass<2; fPass++) { + Int_t &constraint=fConstraint[fPass]; if (constraint<0) continue; + for (Int_t i=0; ifReconstructed) continue; //this track was already "succesfully" reconstructed + if ( (TMath::Abs(t->GetD(GetX(),GetY())) >2.) && fConstraint[fPass]) continue; + if ( (TMath::Abs(t->GetZat(GetX())-GetZ())>2.) && fConstraint[fPass]) continue; + + Int_t tpcLabel=t->GetLabel(); //save the TPC track label + + ResetTrackToFollow(*t); + ResetBestTrack(); + + for (FollowProlongation(); fISetLabel(tpcLabel); + besttrack->CookdEdx(); + besttrack->fFakeRatio=1.; + CookLabel(besttrack,0.); //For comparison only + // besttrack->UpdateESDtrack(AliESDtrack::kITSin); + + if (besttrack->GetChi2()/besttrack->GetNumberOfClusters()>3.5){ + if ( (TMath::Abs(besttrack->GetD(GetX(),GetY()))>0.4) && fConstraint[fPass]) { + CompressTrackHypothesys(fCurrentEsdTrack,0.0,0); + continue; + } + if ( (TMath::Abs(besttrack->GetZat(GetX()) -GetZ() )>0.4) && fConstraint[fPass]){ + CompressTrackHypothesys(fCurrentEsdTrack,0.0,0); + continue; + } + } + + //delete itsTracks.RemoveAt(i); + t->fReconstructed = kTRUE; + ntrk++; + + } + } + + for (Int_t i=0; iGetLabel(); //save the TPC track label + AliITStrackV2 * besttrack = GetBestHypothesysMIP(fCurrentEsdTrack,t); + if (!besttrack) continue; + + besttrack->SetLabel(tpcLabel); + besttrack->CookdEdx(); + besttrack->fFakeRatio=1.; + CookLabel(besttrack,0.); //For comparison only + besttrack->UpdateESDtrack(AliESDtrack::kITSin); } + // - if (!out->IsOpen()) { - Error("Clusters2Tracks","file for ITS tracks is not open !\n"); - return 2; + itsTracks.Delete(); + // + Int_t entries = fTrackHypothesys.GetEntriesFast(); + for (Int_t ientry=0;ientryDelete(); + delete fTrackHypothesys.RemoveAt(ientry); } - in->cd(); - - Char_t tname[100]; + fTrackHypothesys.Delete(); + Info("Clusters2Tracks","Number of prolonged tracks: %d\n",ntrk); + + return 0; +} + +Int_t AliITStrackerV2::Clusters2Tracks(TTree *tpcTree, TTree *itsTree) { + //-------------------------------------------------------------------- + // This functions reconstructs ITS tracks + // The clusters must be already loaded ! + //-------------------------------------------------------------------- Int_t nentr=0; TObjArray itsTracks(15000); + Warning("Clusters2Tracks(TTree *, TTree *)", + "Will be removed soon ! Use Clusters2Tracks(AliESD *) instead."); + {/* Read TPC tracks */ - sprintf(tname,"TreeT_TPC_%d",GetEventNumber()); - TTree *tpcTree=(TTree*)in->Get(tname); - if (!tpcTree) { - Error("Clusters2Tracks","can't get a tree with TPC tracks !\n"); - return 3; - } AliTPCtrack *itrack=new AliTPCtrack; + TBranch *branch=tpcTree->GetBranch("tracks"); + if (!branch) { + Error("Clusters2Tracks","Can't get the branch !"); + return 1; + } tpcTree->SetBranchAddress("tracks",&itrack); nentr=(Int_t)tpcTree->GetEntries(); @@ -214,28 +368,30 @@ Int_t AliITStrackerV2::Clusters2Tracks(const TFile *inp, TFile *out) { delete t; continue; } - if (TMath::Abs(t->GetD())>4) continue; + if (TMath::Abs(t->GetD())>4) { + delete t; + continue; + } if (CorrectForDeadZoneMaterial(t)!=0) { Warning("Clusters2Tracks", "failed to correct for the material in the dead zone !\n"); + delete t; continue; } itsTracks.AddLast(t); } - delete tpcTree; //Thanks to Mariana Bondila delete itrack; } itsTracks.Sort(); nentr=itsTracks.GetEntriesFast(); - out->cd(); - sprintf(tname,"TreeT_ITS_%d",GetEventNumber()); - TTree itsTree(tname,"Tree with ITS tracks"); AliITStrackV2 *otrack=&fBestTrack; - itsTree.Branch("tracks","AliITStrackV2",&otrack,32000,0); + TBranch *branch=itsTree->GetBranch("tracks"); + if (!branch) itsTree->Branch("tracks","AliITStrackV2",&otrack,32000,3); + else branch->SetAddress(&otrack); for (fPass=0; fPass<2; fPass++) { Int_t &constraint=fConstraint[fPass]; if (constraint<0) continue; @@ -251,241 +407,170 @@ Int_t AliITStrackerV2::Clusters2Tracks(const TFile *inp, TFile *out) { while (TakeNextProlongation()) FollowProlongation(); } - if (fBestTrack.GetNumberOfClusters() < kMaxLayer-kLayersToSkip)continue; + if (fBestTrack.GetNumberOfClusters() == 0) continue; if (fConstraint[fPass]) { - if (!RefitAt(3.7, t, &fBestTrack)) continue; + ResetTrackToFollow(*t); + if (!RefitAt(3.7, &fTrackToFollow, &fBestTrack)) continue; + ResetBestTrack(); } fBestTrack.SetLabel(tpcLabel); fBestTrack.CookdEdx(); CookLabel(&fBestTrack,0.); //For comparison only - itsTree.Fill(); - UseClusters(&fBestTrack); + itsTree->Fill(); + //UseClusters(&fBestTrack); delete itsTracks.RemoveAt(i); } } - nentr=(Int_t)itsTree.GetEntries(); + nentr=(Int_t)itsTree->GetEntries(); Info("Clusters2Tracks","Number of prolonged tracks: %d\n",nentr); - itsTree.Write(); - itsTracks.Delete(); - UnloadClusters(); - - savedir->cd(); return 0; } -Int_t AliITStrackerV2::PropagateBack(const TFile *inp, TFile *out) { +Int_t AliITStrackerV2::PropagateBack(AliESD *event) { //-------------------------------------------------------------------- - //This functions propagates reconstructed ITS tracks back + // This functions propagates reconstructed ITS tracks back + // The clusters must be loaded ! //-------------------------------------------------------------------- - TFile *in=(TFile*)inp; - TDirectory *savedir=gDirectory; - - if (LoadClusters()!=0) return 1; - - if (!in->IsOpen()) { - Error("PropagateBack","file with ITS tracks is not open !\n"); - return 1; - } - - if (!out->IsOpen()) { - Error("PropagateBack","file for back propagated ITS tracks is not open !\n"); - return 2; - } - - in->cd(); - - Char_t tname[100]; - sprintf(tname,"TreeT_ITS_%d",GetEventNumber()); - TTree *itsTree=(TTree*)in->Get(tname); - if (!itsTree) { - Error("PropagateBack","can't get a tree with ITS tracks !\n"); - return 3; - } - AliITStrackV2 *itrack=new AliITStrackV2; - itsTree->SetBranchAddress("tracks",&itrack); - - out->cd(); - - sprintf(tname,"TreeT_ITSb_%d",GetEventNumber()); - TTree backTree(tname,"Tree with back propagated ITS tracks"); - AliTPCtrack *otrack=0; - backTree.Branch("tracks","AliTPCtrack",&otrack,32000,2); + Int_t nentr=event->GetNumberOfTracks(); + Info("PropagateBack", "Number of ESD tracks: %d\n", nentr); + + Int_t ntrk=0; + for (Int_t i=0; iGetTrack(i); + + if ((esd->GetStatus()&AliESDtrack::kITSin)==0) continue; + if (esd->GetStatus()&AliESDtrack::kITSout) continue; + + AliITStrackV2 *t=0; + try { + t=new AliITStrackV2(*esd); + } catch (const Char_t *msg) { + Warning("PropagateBack",msg); + delete t; + continue; + } - Int_t nentr=(Int_t)itsTree->GetEntries(); - Int_t i; - for (i=0; iGetEvent(i); - Int_t itsLabel=itrack->GetLabel(); //save the ITS track label - ResetTrackToFollow(*itrack); - - // propagete to vertex [SR, GSI 17.02.2003] - fTrackToFollow.PropagateTo(3.,0.0028,65.19); - fTrackToFollow.PropagateToVertex(); - - // Start Time measurement [SR, GSI 17.02.2003] - fTrackToFollow.StartTimeIntegral(); - fTrackToFollow.PropagateTo(3.,-0.0028,65.19); - // + ResetTrackToFollow(*t); - itrack->ResetCovariance(); itrack->ResetClusters(); - if (!RefitAt(49.,itrack,&fTrackToFollow)) continue; + // propagete to vertex [SR, GSI 17.02.2003] + // Start Time measurement [SR, GSI 17.02.2003], corrected by I.Belikov + if (fTrackToFollow.PropagateTo(3.,0.0028,65.19)) { + if (fTrackToFollow.PropagateToVertex()) { + fTrackToFollow.StartTimeIntegral(); + } + fTrackToFollow.PropagateTo(3.,-0.0028,65.19); + } - if (CorrectForDeadZoneMaterial(&fTrackToFollow)!=0) { - Warning("PropagateBack", - "failed to correct for the material in the dead zone !\n"); - continue; - } - - fTrackToFollow.SetLabel(itsLabel); - otrack=new AliTPCtrack(fTrackToFollow,fTrackToFollow.GetAlpha()); - backTree.Fill(); delete otrack; - UseClusters(&fTrackToFollow); + fTrackToFollow.ResetCovariance(); fTrackToFollow.ResetClusters(); + if (RefitAt(49.,&fTrackToFollow,t)) { + if (CorrectForDeadZoneMaterial(&fTrackToFollow)!=0) { + Warning("PropagateBack", + "failed to correct for the material in the dead zone !\n"); + delete t; + continue; + } + fTrackToFollow.SetLabel(t->GetLabel()); + //fTrackToFollow.CookdEdx(); + CookLabel(&fTrackToFollow,0.); //For comparison only + fTrackToFollow.UpdateESDtrack(AliESDtrack::kITSout); + //UseClusters(&fTrackToFollow); + ntrk++; + } + delete t; } - i=(Int_t)backTree.GetEntries(); - backTree.Write(); - - Info("PropagateBack","Number of ITS tracks: %d\n",nentr); - Info("PropagateBack","Number of back propagated ITS tracks: %d\n",i); - - delete itrack; - delete itsTree; //Thanks to Mariana Bondila - UnloadClusters(); - - savedir->cd(); + Info("PropagateBack","Number of back propagated ITS tracks: %d\n",ntrk); return 0; } -Int_t AliITStrackerV2::RefitInward(const TFile *inp, TFile *out) { +Int_t AliITStrackerV2::RefitInward(AliESD *event) { //-------------------------------------------------------------------- // This functions refits ITS tracks using the // "inward propagated" TPC tracks + // The clusters must be loaded ! //-------------------------------------------------------------------- - TFile *in=(TFile*)inp; - TDirectory *savedir=gDirectory; + Int_t nentr=event->GetNumberOfTracks(); + Info("RefitInward", "Number of ESD tracks: %d\n", nentr); - if (LoadClusters()!=0) return 1; + Int_t ntrk=0; + for (Int_t i=0; iGetTrack(i); - if (!in->IsOpen()) { - Error("RefitInward","file with inward TPC tracks is not open !\n"); - return 2; - } - - if (!out->IsOpen()) { - Error("RefitInward","file for inward ITS tracks is not open !\n"); - return 3; - } - - Int_t i; + if ((esd->GetStatus()&AliESDtrack::kITSout) == 0) continue; + if (esd->GetStatus()&AliESDtrack::kITSrefit) continue; + if (esd->GetStatus()&AliESDtrack::kTPCout) + if ((esd->GetStatus()&AliESDtrack::kTPCrefit)==0) continue; - //LUT used for the track matching (S.Radomski's idea) - const Int_t nLab = 400000; // limit on the number of track labels - Int_t lut[nLab]; - for(i=0; iGet(tname); - if (!itsTree) { - Error("RefitInward","can't get a tree with ITS tracks !\n"); - return 3; - } - AliITStrackV2 *itrack=new AliITStrackV2; - itsTree->SetBranchAddress("tracks",&itrack); - Int_t nits=(Int_t)itsTree->GetEntries(); - - Info("RefitInward","Number of ITS tracks: %d\n",nits); - - for (Int_t i=0; iGetEvent(i); - Int_t lab=TMath::Abs(itrack->GetLabel()); - if (lab < nLab) { - if (lut[lab]>=0) Warning("RefitInward","double track %d\n",lab); - lut[lab]=i; - } else { - Warning("RefitInward","Too big ITS track label: %d\n!",lab); - continue; - } - itsTracks.AddLast(new AliITStrackV2(*itrack)); - } - delete itsTree; - delete itrack; - } - - out->cd(); - - //Create the output tree - sprintf(tname,"TreeT_ITSinward_%d",GetEventNumber()); - TTree outTree(tname,"Tree with inward refitted ITS tracks"); - AliITStrackV2 *otrack=0; - outTree.Branch("tracks","AliITStrackV2",&otrack,32000,0); - - //Get the input tree - sprintf(tname,"tracksTPC_%d",GetEventNumber()); - TTree *tpcTree=(TTree*)in->Get(tname); - if (!tpcTree) { - Error("RefitInward","can't get a tree with TPC tracks !\n"); - return 3; - } - AliTPCtrack *itrack=new AliTPCtrack; - tpcTree->SetBranchAddress("tracks",&itrack); - Int_t ntpc=(Int_t)tpcTree->GetEntries(); - - Info("RefitInward","Number of TPC tracks: %d\n",ntpc); - - for (i=0; iGetEvent(i); AliITStrackV2 *t=0; try { - t=new AliITStrackV2(*itrack); + t=new AliITStrackV2(*esd); } catch (const Char_t *msg) { - Warning("RefitInward",msg); - delete t; - continue; - } - //check if this track was reconstructed in the ITS - Int_t lab=TMath::Abs(t->GetLabel()); - if (lab >= nLab) { - Warning("RefitInward","Too big TPC track label: %d\n!",lab); - continue; + Warning("RefitInward",msg); + delete t; + continue; } - Int_t idx=lut[lab]; - if (idx<0) continue; //no prolongation in the ITS for this track - + if (CorrectForDeadZoneMaterial(t)!=0) { Warning("RefitInward", "failed to correct for the material in the dead zone !\n"); + delete t; continue; } + ResetTrackToFollow(*t); + fTrackToFollow.ResetClusters(); + + if ((esd->GetStatus()&AliESDtrack::kTPCin)==0) + fTrackToFollow.ResetCovariance(); + //Refitting... - otrack=(AliITStrackV2*)itsTracks.UncheckedAt(idx); - if (!RefitAt(3.7, t, otrack)) continue; - otrack->SetLabel(itrack->GetLabel()); //For comparison only - otrack->CookdEdx(); - CookLabel(otrack,0.); //For comparison only - outTree.Fill(); + if (RefitAt(3.7, &fTrackToFollow, t)) { + fTrackToFollow.SetLabel(t->GetLabel()); + fTrackToFollow.CookdEdx(); + CookLabel(&fTrackToFollow,0.0); //For comparison only + + if (fTrackToFollow.PropagateTo(3.,0.0028,65.19)) {//The beam pipe + Double_t a=fTrackToFollow.GetAlpha(); + Double_t cs=TMath::Cos(a),sn=TMath::Sin(a); + Double_t xv= GetX()*cs + GetY()*sn; + Double_t yv=-GetX()*sn + GetY()*cs; + + Double_t c=fTrackToFollow.GetC(), snp=fTrackToFollow.GetSnp(); + Double_t x=fTrackToFollow.GetX(), y=fTrackToFollow.GetY(); + Double_t tgfv=-(c*(x-xv)-snp)/(c*(y-yv) + TMath::Sqrt(1.-snp*snp)); + Double_t fv=TMath::ATan(tgfv); + + cs=TMath::Cos(fv); sn=TMath::Sin(fv); + x = xv*cs + yv*sn; + yv=-xv*sn + yv*cs; xv=x; + + if (fTrackToFollow.Propagate(fv+a,xv)) { + fTrackToFollow.UpdateESDtrack(AliESDtrack::kITSrefit); + //UseClusters(&fTrackToFollow); + { + AliITSclusterV2 c; c.SetY(yv); c.SetZ(GetZ()); + c.SetSigmaY2(GetSigmaY()*GetSigmaY()); + c.SetSigmaZ2(GetSigmaZ()*GetSigmaZ()); + Double_t chi2=fTrackToFollow.GetPredictedChi2(&c); + if (chi2cd(); + Info("RefitInward","Number of refitted tracks: %d\n",ntrk); return 0; } @@ -496,7 +581,7 @@ AliCluster *AliITStrackerV2::GetCluster(Int_t index) const { //-------------------------------------------------------------------- Int_t l=(index & 0xf0000000) >> 28; Int_t c=(index & 0x0fffffff) >> 00; - return fLayers[l].GetCluster(c); + return fgLayers[l].GetCluster(c); } @@ -504,23 +589,21 @@ void AliITStrackerV2::FollowProlongation() { //-------------------------------------------------------------------- //This function finds a track prolongation //-------------------------------------------------------------------- - Int_t tryAgain=kLayersToSkip; - - while (fI) { + while (fI>fLastLayerToTrackTo) { Int_t i=fI-1; - AliITSlayer &layer=fLayers[i]; + AliITSlayer &layer=fgLayers[i]; AliITStrackV2 &track=fTracks[i]; Double_t r=layer.GetR(); if (i==3 || i==1) { - Double_t rs=0.5*(fLayers[i+1].GetR() + r); + Double_t rs=0.5*(fgLayers[i+1].GetR() + r); Double_t d=0.0034, x0=38.6; if (i==1) {rs=9.; d=0.0097; x0=42;} if (!fTrackToFollow.PropagateTo(rs,d,x0)) { //Warning("FollowProlongation","propagation failed !\n"); - break; + return; } } @@ -528,13 +611,14 @@ void AliITStrackerV2::FollowProlongation() { Double_t x,y,z; if (!fTrackToFollow.GetGlobalXYZat(r,x,y,z)) { //Warning("FollowProlongation","failed to estimate track !\n"); - break; + return; } Double_t phi=TMath::ATan2(y,x); + Int_t idet=layer.FindDetectorIndex(phi,z); if (idet<0) { //Warning("FollowProlongation","failed to find a detector !\n"); - break; + return; } //propagate to the intersection @@ -542,7 +626,7 @@ void AliITStrackerV2::FollowProlongation() { phi=det.GetPhi(); if (!fTrackToFollow.Propagate(phi,det.GetR())) { //Warning("FollowProlongation","propagation failed !\n"); - break; + return; } fTrackToFollow.SetDetectorIndex(idet); @@ -560,16 +644,20 @@ void AliITStrackerV2::FollowProlongation() { if (dz < 0.5*TMath::Abs(track.GetTgl())) dz=0.5*TMath::Abs(track.GetTgl()); if (dz > kMaxRoad) { //Warning("FollowProlongation","too broad road in Z !\n"); - break; + return; } - if (TMath::Abs(fTrackToFollow.GetZ()-GetZ()) > r+dz) break; + // if (TMath::Abs(fTrackToFollow.GetZ()-GetZ()) > r+dz) return; //Double_t dy=4*TMath::Sqrt(track.GetSigmaY2() + kSigmaY2[i]); + if (TMath::Abs(track.GetSnp()>kMaxSnp)) { + fI--; + continue; // MI + } if (dy < 0.5*TMath::Abs(track.GetSnp())) dy=0.5*TMath::Abs(track.GetSnp()); if (dy > kMaxRoad) { //Warning("FollowProlongation","too broad road in Y !\n"); - break; + return; } Double_t zmin=track.GetZ() - dz; @@ -580,24 +668,49 @@ void AliITStrackerV2::FollowProlongation() { fI--; //take another prolongation - if (!TakeNextProlongation()) if (!tryAgain--) break; - tryAgain=kLayersToSkip; - + if (!TakeNextProlongation()){ + //skipped++; + fTrackToFollow.fNSkipped++; + if (fLayersNotToSkip[fI]||fTrackToFollow.fNSkipped>1) return; + } + if (fTrackToFollow.fNUsed>1) return; + if (fTrackToFollow.fNUsed+fTrackToFollow.fNSkipped>1) return; + if ( (fI<3) && ( fTrackToFollow.GetChi2()/fTrackToFollow.GetNumberOfClusters()>kChi2PerCluster)) return; } //deal with the best track Int_t ncl=fTrackToFollow.GetNumberOfClusters(); Int_t nclb=fBestTrack.GetNumberOfClusters(); - if (ncl) - if (ncl >= nclb) { - Double_t chi2=fTrackToFollow.GetChi2(); - if (chi2/ncl < kChi2PerCluster) { - if (ncl > nclb || chi2 < fBestTrack.GetChi2()) { - ResetBestTrack(); + if (ncl){ + if (ncl<4) return; + if ( (ncl<6) && (fTrackToFollow.GetChi2()/float(ncl))>3) return; + if (fTrackToFollow.GetChi2()/ncl>5.5) return; + fTrackToFollow.CookdEdx(); + if (fTrackToFollow.fESDtrack->GetTPCsignal()>80.) + if ((fTrackToFollow.GetdEdx()/fTrackToFollow.fESDtrack->GetTPCsignal())<0.35){ + // mismatch in dedx + return; + } + // + fTrackToFollow.SetLabel(fTrackToFollow.fESDtrack->GetLabel()); + CookLabel(&fTrackToFollow,0.); // + // + if ( (nclb>3) && ((fTrackToFollow.GetChi2()/ncl)<(3*fBestTrack.GetChi2()/(nclb)))) + AddTrackHypothesys(new AliITStrackV2(fTrackToFollow), fCurrentEsdTrack); + else + if (ncl>3) AddTrackHypothesys(new AliITStrackV2(fTrackToFollow), fCurrentEsdTrack); + if (ncl >= nclb) { + Double_t chi2=fTrackToFollow.GetChi2(); + if (chi2/ncl < kChi2PerCluster) { + if (ncl > nclb ) { + ResetBestTrack(); } - } + if ( (ncl == nclb) && chi2 < fBestTrack.GetChi2()) { + ResetBestTrack(); + } + } + } } - } Int_t AliITStrackerV2::TakeNextProlongation() { @@ -606,7 +719,7 @@ Int_t AliITStrackerV2::TakeNextProlongation() { // // dEdx analysis by: Boris Batyunya, JINR, Boris.Batiounia@cern.ch //-------------------------------------------------------------------- - AliITSlayer &layer=fLayers[fI]; + AliITSlayer &layer=fgLayers[fI]; ResetTrackToFollow(fTracks[fI]); Double_t dz=7*TMath::Sqrt(fTrackToFollow.GetSigmaZ2() + kSigmaZ2[fI]); @@ -621,7 +734,6 @@ Int_t AliITStrackerV2::TakeNextProlongation() { Double_t chi2=12345.; while ((c=layer.GetNextCluster(ci))!=0) { Int_t idet=c->GetDetectorIndex(); - if (fTrackToFollow.GetDetectorIndex()!=idet) { const AliITSdetector &det=layer.GetDetector(idet); ResetTrackToFollow(fTracks[fI]); @@ -641,12 +753,23 @@ Int_t AliITStrackerV2::TakeNextProlongation() { if (chi2>=kMaxChi2) return 0; if (!c) return 0; + if (c->IsUsed()&&c->GetNy()<5) { //shared factor + chi2+=1; + chi2*=2*(1./(TMath::Max(c->GetNy(),1))); + } + if (c->GetQ()==0){ //dead zone virtual cluster + chi2*=4.; + chi2+=20; + fTrackToFollow.fNUsed++; + return 1; + } + //if ((fI<2)&&chi2>kMaxChi2In) return 0; if (!fTrackToFollow.Update(c,chi2,(fI<<28)+ci)) { //Warning("TakeNextProlongation","filtering failed !\n"); return 0; } - + if (c->IsUsed()&&c->GetNy()<5) fTrackToFollow.fNUsed++; if (fTrackToFollow.GetNumberOfClusters()>1) if (TMath::Abs(fTrackToFollow.GetD())>4) return 0; @@ -663,6 +786,15 @@ Int_t AliITStrackerV2::TakeNextProlongation() { Double_t d=GetEffectiveThickness(0,0); //Think of this !!!! Double_t xyz[]={GetX(),GetY(),GetZ()}; Double_t ers[]={GetSigmaX(),GetSigmaY(),GetSigmaZ()}; + Double_t deltad = TMath::Abs(fTrackToFollow.GetD(GetX(),GetY())); + Double_t deltaz = TMath::Abs(fTrackToFollow.GetZat(GetX())-GetZ()); + + if ( (fI==4) && (deltad>2.0 || deltaz>1.5)) return 0; // don't "improve" secondaries + if ( (fI==3) && (deltad>1.5 || deltaz>0.9)) return 0; // don't "improve" secondaries + if ( (fI==2) && (deltad>0.9 || deltaz>0.6)) return 1; // don't "improve" secondaries + if ( (fI==1) && (deltad>0.3 || deltaz>0.3)) return 1; // don't "improve" secondaries + if ( (fI==0) && (deltad>0.1 || deltaz>0.1)) return 1; // don't "improve" secondaries + fTrackToFollow.Improve(d,xyz,ers); } @@ -699,6 +831,7 @@ AliITStrackerV2::AliITSlayer::~AliITSlayer() { //-------------------------------------------------------------------- delete[] fDetectors; for (Int_t i=0; iGetZ() > fZmax) break; - if (c->IsUsed()) continue; + // if (c->IsUsed()) continue; const AliITSdetector &det=GetDetector(c->GetDetectorIndex()); Double_t y=fR*det.GetPhi() + c->GetY(); @@ -794,7 +935,7 @@ FindDetectorIndex(Double_t phi, Double_t z) const { //-------------------------------------------------------------------- //This function finds the detector crossed by the track //-------------------------------------------------------------------- - Double_t dphi=phi-fPhiOffset; + Double_t dphi=-(phi-fPhiOffset); if (dphi < 0) dphi += 2*TMath::Pi(); else if (dphi >= 2*TMath::Pi()) dphi -= 2*TMath::Pi(); Int_t np=Int_t(dphi*fNladders*0.5/TMath::Pi()+0.5); @@ -919,10 +1060,10 @@ Double_t AliITStrackerV2::GetEffectiveThickness(Double_t y,Double_t z) const Double_t d=0.0028*3*3; //beam pipe Double_t x0=0; - Double_t xn=fLayers[fI].GetR(); + Double_t xn=fgLayers[fI].GetR(); for (Int_t i=0; i1) { @@ -931,7 +1072,7 @@ Double_t AliITStrackerV2::GetEffectiveThickness(Double_t y,Double_t z) const } if (fI>3) { - Double_t xi=0.5*(fLayers[3].GetR()+fLayers[4].GetR()); + Double_t xi=0.5*(fgLayers[3].GetR()+fgLayers[4].GetR()); d+=0.0034*xi*xi; } @@ -961,23 +1102,22 @@ Int_t AliITStrackerV2::AliITSlayer::InRoad() const { } Bool_t -AliITStrackerV2::RefitAt(Double_t x,const AliITStrackV2 *s,AliITStrackV2 *ot) { +AliITStrackerV2::RefitAt(Double_t xx,AliITStrackV2 *t,const AliITStrackV2 *c) { //-------------------------------------------------------------------- - // This function refits a track at a given position + // This function refits the track "t" at the position "x" using + // the clusters from "c" //-------------------------------------------------------------------- - AliITStrackV2 save(*ot), *t=&save; Int_t index[kMaxLayer]; Int_t k; for (k=0; kGetNumberOfClusters(); + Int_t nc=c->GetNumberOfClusters(); for (k=0; kGetClusterIndex(k),nl=(idx&0xf0000000)>>28; + Int_t idx=c->GetClusterIndex(k),nl=(idx&0xf0000000)>>28; index[nl]=idx; } - t->~AliITStrackV2(); new (t) AliITStrackV2(*s); Int_t from, to, step; - if (x > t->GetX()) { + if (xx > t->GetX()) { from=0; to=kMaxLayer; step=+1; } else { @@ -986,13 +1126,13 @@ AliITStrackerV2::RefitAt(Double_t x,const AliITStrackV2 *s,AliITStrackV2 *ot) { } for (Int_t i=from; i != to; i += step) { - AliITSlayer &layer=fLayers[i]; + AliITSlayer &layer=fgLayers[i]; Double_t r=layer.GetR(); { Double_t hI=i-0.5*step; if (TMath::Abs(hI-1.5)<0.01 || TMath::Abs(hI-3.5)<0.01) { - Double_t rs=0.5*(fLayers[i-step].GetR() + r); + Double_t rs=0.5*(fgLayers[i-step].GetR() + r); Double_t d=0.0034, x0=38.6; if (TMath::Abs(hI-1.5)<0.01) {rs=9.; d=0.0097; x0=42;} if (!t->PropagateTo(rs,-step*d,x0)) { @@ -1039,8 +1179,12 @@ AliITStrackerV2::RefitAt(Double_t x,const AliITStrackV2 *s,AliITStrackV2 *ot) { t->SetDetectorIndex(idet); } Double_t chi2=t->GetPredictedChi2(c); - if (chi2PropagateTo(x,0.,0.)) return kFALSE; - ot->~AliITStrackV2(); new (ot) AliITStrackV2(*t); + if (!t->PropagateTo(xx,0.,0.)) return kFALSE; return kTRUE; } + +Float_t *AliITStrackerV2::GetWeight(Int_t index) { + //-------------------------------------------------------------------- + // Return pointer to a given cluster + //-------------------------------------------------------------------- + Int_t l=(index & 0xf0000000) >> 28; + Int_t c=(index & 0x0fffffff) >> 00; + return fgLayers[l].GetWeight(c); +} + + + void AliITStrackerV2::UseClusters(const AliKalmanTrack *t, Int_t from) const { //-------------------------------------------------------------------- // This function marks clusters assigned to the track @@ -1105,3 +1260,537 @@ void AliITStrackerV2::UseClusters(const AliKalmanTrack *t, Int_t from) const { if (c->GetSigmaZ2()>0.1) c->Use(); } + + +void AliITStrackerV2::AddTrackHypothesys(AliITStrackV2 * track, Int_t esdindex) +{ + //------------------------------------------------------------------ + // add track to the list of hypothesys + //------------------------------------------------------------------ + + if (esdindex>=fTrackHypothesys.GetEntriesFast()) fTrackHypothesys.Expand(esdindex*2+10); + // + TObjArray * array = (TObjArray*) fTrackHypothesys.At(esdindex); + if (!array) { + array = new TObjArray(10); + fTrackHypothesys.AddAt(array,esdindex); + } + array->AddLast(track); +} + +void AliITStrackerV2::SortTrackHypothesys(Int_t esdindex, Float_t likelihoodlevel) +{ + //------------------------------------------------------------------- + // compress array of track hypothesys + // keep only maxsize best hypothesys + //------------------------------------------------------------------- + if (esdindex>fTrackHypothesys.GetEntriesFast()) return; + if (! (fTrackHypothesys.At(esdindex)) ) return; + TObjArray * array = (TObjArray*) fTrackHypothesys.At(esdindex); + Int_t entries = array->GetEntriesFast(); + + Float_t * chi2 = new Float_t[entries]; + Float_t * probability = new Float_t[entries]; + Float_t sumprobabilityall=0; + Int_t * index = new Int_t[entries]; + // + // + for (Int_t itrack=0;itrackGetEntriesFast();itrack++){ + AliITStrackV2 * track = (AliITStrackV2*)array->At(itrack); + // + if (track && track->GetNumberOfClusters()>(track->fNUsed+track->fNSkipped)){ + // + chi2[itrack] = track->GetChi2()/(track->GetNumberOfClusters()-track->fNUsed-track->fNSkipped); + if (track->fESDtrack) + if (track->fESDtrack->GetTPCsignal()>80){ + track->CookdEdx(); + if ((track->GetdEdx()/track->fESDtrack->GetTPCsignal())<0.4){ + Float_t factor = 2.+10.*(0.6-track->GetdEdx()/track->fESDtrack->GetTPCsignal()); + chi2[itrack]*= factor; //mismatch in dEdx + } + } + } + else + chi2[itrack] = 10000000.; + probability[itrack] = 1./(0.3+chi2[itrack]); + sumprobabilityall+=probability[itrack]; + } + + TMath::Sort(entries,chi2,index,kFALSE); + TObjArray * newarray = new TObjArray(); + Float_t sumprobability = 0.; + for (Int_t i=0;iAt(index[i]); + if (!track) break; + if (chi2[index[i]]<30){ + newarray->AddLast(array->RemoveAt(index[i])); + track->fChi2MIP[0] = chi2[index[i]]; // base chi 2 + sumprobability+= probability[index[i]]/sumprobabilityall; + if (sumprobability> likelihoodlevel) break; + } + else{ + delete array->RemoveAt(index[i]); + } + } + + array->Delete(); + delete fTrackHypothesys.RemoveAt(esdindex); + fTrackHypothesys.AddAt(newarray,esdindex); + + delete [] chi2; + delete [] index; + +} + + +void AliITStrackerV2::CompressTrackHypothesys(Int_t esdindex, Float_t likelihoodlevel, Int_t maxsize) +{ + // + // + if (esdindex>fTrackHypothesys.GetEntriesFast()) return; + if (! (fTrackHypothesys.At(esdindex)) ) return; + TObjArray * array = (TObjArray*) fTrackHypothesys.At(esdindex); + Int_t entries = array->GetEntriesFast(); + // + if (likelihoodlevel>0.000001){ + Float_t *probability = new Float_t[entries]; + for (Int_t i=0;iAt(itrack); + probability[itrack]=0; + if (!track) continue; + probability[itrack] = 1./(0.3+track->fChi2MIP[0]); + sumprobabilityall += probability[itrack]; + // + } + if (sumprobabilityall<=0.000000000000001){ + return; + } + for (Int_t itrack=0;itrackAt(index[i]); + if (!track) continue; + newarray->AddLast(array->RemoveAt(index[i])); + sumprobability+= probability[index[i]]; + if (sumprobability> likelihoodlevel) break; + if (i>maxsize) break; + } + + array->Delete(); + delete fTrackHypothesys.RemoveAt(esdindex); + fTrackHypothesys.AddAt(newarray,esdindex); + // + delete []index; + delete []probability; + } + else{ + array->Delete(); + delete fTrackHypothesys.RemoveAt(esdindex); + } +} + + +AliITStrackV2 * AliITStrackerV2::GetBestHypothesys(Int_t esdindex, AliITStrackV2 * original, Int_t checkmax) +{ + //------------------------------------------------------------- + // try to find best hypothesy + // currently - minimal chi2 of track+backpropagated track+matching to the tpc track + //------------------------------------------------------------- + if (fTrackHypothesys.GetEntriesFast()<=esdindex) return 0; + TObjArray * array = (TObjArray*) fTrackHypothesys.At(esdindex); + if (!array) return 0; + Int_t entries = array->GetEntriesFast(); + if (!entries) return 0; + Float_t minchi2 = 100000; + Int_t maxn = 3; + AliITStrackV2 * besttrack=0; + Int_t accepted =0; + Int_t maxindex=0; + // + Float_t sumz2=0; + Float_t sumy2=0; + Float_t sumall=0; + for (Int_t itrack=0;itrackAt(itrack); + if (!track) continue; + sumall++; + sumz2+=track->GetSigmaZ2(); + sumy2+=track->GetSigmaY2(); + } + sumz2/=sumall; + sumy2/=sumall; + + Float_t dedxmismatch=1; + for (Int_t i=0;iAt(i); + if (!track) continue; + track->fChi2MIP[1] = 1000000; + track->fChi2MIP[2] = 1000000; + + if ( (track->GetNumberOfClusters()-track->fNSkipped-track->fNUsed)<2) continue; + // + if (track->fESDtrack) + if (track->fESDtrack->GetTPCsignal()>80){ + track->CookdEdx(); + if ((track->GetdEdx()/track->fESDtrack->GetTPCsignal())<0.4){ + //mismatch in dEdx + dedxmismatch= 2.+ 10.*(0.6-track->GetdEdx()/track->fESDtrack->GetTPCsignal()); + } + } + + // track->SetLabel(original->GetLabel()); + //CookLabel(track,0.0); + //if (track->GetFakeRatio()>0.01) continue; + // + // + // backtrack + AliITStrackV2 * backtrack = new AliITStrackV2(*track); + backtrack->ResetCovariance(); + backtrack->ResetClusters(); + Double_t x = original->GetX(); + if (!RefitAt(x,backtrack,track)){ + delete backtrack; + delete array->RemoveAt(i); + continue; + } + if ( (backtrack->GetChi2() / float(backtrack->GetNumberOfClusters()-track->fNSkipped-track->fNUsed-0.5))>6) + { + delete backtrack; + delete array->RemoveAt(i); + continue; + } + Double_t deltac = backtrack->GetC()-original->GetC(); + Double_t deltatgl = backtrack->GetTgl()-original->GetTgl(); + // + Double_t poolc2 = (deltac*deltac)/(original->fC44+backtrack->fC44); + Double_t pooltgl2 = (deltatgl*deltatgl)/(original->fC33+backtrack->fC33); + if ((poolc2+pooltgl2)>32){ //4 sigma + delete backtrack; + delete array->RemoveAt(i); + continue; + } + //Double_t bpoolc = (deltac*deltac)/(original->fC44); + //Double_t bpooltgl = (deltatgl*deltatgl)/(original->fC33); + + // + //forward track - without constraint + AliITStrackV2 * forwardtrack = new AliITStrackV2(*original); + // forwardtrack->ResetCovariance(); + forwardtrack->ResetClusters(); + x = track->GetX(); + if (!RefitAt(x,forwardtrack,track)){ + delete forwardtrack; + delete backtrack; + delete array->RemoveAt(i); + continue; + } + if ( (forwardtrack->GetChi2()/float(forwardtrack->GetNumberOfClusters()-track->fNSkipped-track->fNUsed))>6) + { + delete forwardtrack; + delete array->RemoveAt(i); + continue; + } + // + accepted++; + if (accepted>checkmax){ + delete backtrack; + delete forwardtrack; + break; + } + Double_t chi2 = (backtrack->GetChi2()/(backtrack->GetNumberOfClusters()-1-track->fNSkipped-track->fNUsed)+ + forwardtrack->GetChi2()/(forwardtrack->GetNumberOfClusters()-track->fNSkipped-track->fNUsed)); + // bpoolc+bpooltgl; + // chi2 *= (forwardtrack->GetSigmaZ2()/sumz2+forwardtrack->GetSigmaY2()/sumy2); + chi2 *= dedxmismatch; + // + // + track->fChi2MIP[1] = backtrack->GetChi2()/(backtrack->GetNumberOfClusters()-1-track->fNSkipped-track->fNUsed); + track->fChi2MIP[2] = forwardtrack->GetChi2()/(forwardtrack->GetNumberOfClusters()-track->fNSkipped-track->fNUsed); + track->fChi2MIP[3] = poolc2+pooltgl2; + // + + if (track->GetNumberOfClusters()>maxn){ + besttrack = new AliITStrackV2(*forwardtrack); + maxn = track->GetNumberOfClusters(); + minchi2 = chi2; + delete backtrack; + delete forwardtrack; + continue; + } + // + if (chi2 < minchi2){ + besttrack = new AliITStrackV2(*forwardtrack); + minchi2 = chi2; + } + delete backtrack; + delete forwardtrack; + } + // + // + if (!besttrack || besttrack->GetNumberOfClusters()<4) { + return 0; + } + + // + besttrack->SetLabel(original->GetLabel()); + CookLabel(besttrack,0.0); + // + // calculate "weight of the cluster" + // + + { + //sign usage information for clusters + Int_t clusterindex[6][100]; + Double_t clusterweight[6][100]; + for (Int_t ilayer=0;ilayer<6;ilayer++) + for (Int_t icluster=0;icluster<100;icluster++){ + clusterindex[ilayer][icluster] = -1; + clusterweight[ilayer][icluster] = 0; + } + //printf("%d\t%d\n",esdindex, entries); + // + Float_t sumchi2=0; + for (Int_t itrack=0;itrackAt(itrack); + if (!track) continue; + if (track->fChi2MIP[1]>1000) continue; //not accepted + sumchi2 +=1./(0.3+track->fChi2MIP[1]+track->fChi2MIP[2]); + } + for (Int_t itrack=0;itrackAt(itrack); + if (!track) continue; + if (track->fChi2MIP[1]>1000) continue; //not accepted + for (Int_t icluster=0;iclusterGetNumberOfClusters();icluster++){ + Int_t tindex = track->GetClusterIndex(icluster); + Int_t ilayer = (tindex & 0xf0000000) >> 28; + if (tindex<0) continue; + Int_t cindex =0; + // + for (cindex=0;cindex<100;cindex++){ + if (clusterindex[ilayer][cindex]<0) break; + if (clusterindex[ilayer][cindex]==tindex) break; + } + if (cindex>100) break; + if (clusterindex[ilayer][cindex]!=tindex) clusterindex[ilayer][cindex] = tindex; + clusterweight[ilayer][cindex]+= (1./(0.3+track->fChi2MIP[1]+track->fChi2MIP[2]))* (1./sumchi2); + Float_t *weight = GetWeight(tindex); + + if (weight){ + *weight+= (1./(0.3+track->fChi2MIP[1]+track->fChi2MIP[2]))* (1./sumchi2); + } + } + } + + if (besttrack->GetChi2()/besttrack->GetNumberOfClusters()>3.5) return besttrack; //don't sign clusters + Int_t current=0; + Double_t deltad = besttrack->GetD(GetX(),GetY()); + Double_t deltaz = besttrack->GetZat(GetX()) - GetZ(); + Double_t deltaprim = TMath::Sqrt(deltad*deltad+deltaz*deltaz); + + for (Int_t icluster=0;iclusterGetNumberOfClusters();icluster++){ + Int_t index = besttrack->GetClusterIndex(icluster); + Int_t ilayer = (index & 0xf0000000) >> 28; + AliITSclusterV2 *c = (AliITSclusterV2*)GetCluster(index); + if (!c) continue; + // + for (Int_t icluster=0;icluster<100;icluster++){ + // sign non "doubt" clusters as used + if (clusterindex[ilayer][icluster]!=index) continue; + // Float_t * weight = GetWeight(index); + //if (weight) if (*weight>1){ + // if (c->IsUsed()) continue; + // c->Use(); + //} + if ( (ilayer*0.2+0.2)GetNy()>4) continue; // don sign cluster + if ( (ilayer>1&&clusterweight[ilayer][icluster]>0.7) || (ilayer<2&&clusterweight[ilayer][icluster]>0.8) ){ + current++; + if (c->IsUsed()) continue; + c->Use(); + } + } + } + } + + // + return besttrack; +} + + +AliITStrackV2 * AliITStrackerV2::GetBestHypothesysMIP(Int_t esdindex, AliITStrackV2 * original) +{ + //------------------------------------------------------------- + // try to find best hypothesy + // currently - minimal chi2 of track+backpropagated track+matching to the tpc track + //------------------------------------------------------------- + if (fTrackHypothesys.GetEntriesFast()<=esdindex) return 0; + TObjArray * array = (TObjArray*) fTrackHypothesys.At(esdindex); + if (!array) return 0; + Int_t entries = array->GetEntriesFast(); + if (!entries) return 0; + AliITStrackV2 * besttrack=0; + // + //sign usage information for clusters + Int_t clusterindex[6][100]; + Double_t clusterweight[6][100]; + for (Int_t ilayer=0;ilayer<6;ilayer++) + for (Int_t icluster=0;icluster<100;icluster++){ + clusterindex[ilayer][icluster] = -1; + clusterweight[ilayer][icluster] = 0; + } + //printf("%d\t%d\n",esdindex, entries); + // + Float_t sumchi2=0; + for (Int_t itrack=0;itrackAt(itrack); + if (!track) continue; + if (track->fChi2MIP[1]>1000) continue; //not accepted - before + sumchi2 +=1./(0.3+track->fChi2MIP[1]+track->fChi2MIP[2]); + } + // + // get cluster weight + for (Int_t itrack=0;itrackAt(itrack); + if (!track) continue; + if (track->fChi2MIP[1]>1000) continue; // track not accepted in previous iterration + for (Int_t icluster=0;iclusterGetNumberOfClusters();icluster++){ + Int_t tindex = track->GetClusterIndex(icluster); + Int_t ilayer = (tindex & 0xf0000000) >> 28; + if (tindex<0) continue; + Int_t cindex =0; + // + for (cindex=0;cindex<100;cindex++){ + if (clusterindex[ilayer][cindex]<0) break; + if (clusterindex[ilayer][cindex]==tindex) break; + } + if (cindex>100) break; + if (clusterindex[ilayer][cindex]!=tindex) clusterindex[ilayer][cindex] = tindex; + clusterweight[ilayer][cindex]+= (1./(0.3+track->fChi2MIP[1]+track->fChi2MIP[2]))* (1./sumchi2); + } + } + // + // get cluster relative sharing - factor + // + // + for (Int_t ilayer=0;ilayer<6;ilayer++) + for (Int_t cindex=0;cindex<100;cindex++){ + if (clusterindex[ilayer][cindex]<0) continue; + Int_t tindex = clusterindex[ilayer][cindex]; + Float_t *weight = GetWeight(tindex); + if (!weight){ + printf("Problem 1\n"); // not existing track + continue; + } + if (*weight<(clusterweight[ilayer][cindex]-0.00001)){ + printf("Problem 2\n"); // not normalized probability + continue; + } + AliITSclusterV2 *c = (AliITSclusterV2*)GetCluster(tindex); + if (c->GetNy()<5){ + clusterweight[ilayer][cindex]/= *weight; + } + else{ + Float_t weight2 = TMath::Max(*weight-0.5*clusterweight[ilayer][cindex],0.0000001); + clusterweight[ilayer][cindex]/= weight2; + if ( clusterweight[ilayer][cindex]>1) clusterweight[ilayer][cindex] =1.; + } + } + // + //take to the account sharing factor + // + Float_t chi2 =10000000; + Float_t sharefactor=0; + Float_t minchi2 = 100000000; + Float_t secchi2 = 100000000; + Float_t norm=0; + for (Int_t itrack=0;itrackAt(itrack); + if (!track) continue; + if (track->fChi2MIP[1]>1000) continue; //not accepted - before + chi2 = track->fChi2MIP[1]; + Float_t newchi2=0; + sharefactor=0; + norm =0; + // + for (Int_t icluster=0;iclusterGetNumberOfClusters();icluster++){ + Int_t tindex = track->GetClusterIndex(icluster); + Int_t ilayer = (tindex & 0xf0000000) >> 28; + if (tindex<0) continue; + Int_t cindex =0; + Float_t cchi2 = (track->fDy[ilayer]*track->fDy[ilayer])/(track->fSigmaY[ilayer]*track->fSigmaY[ilayer]) + + (track->fDz[ilayer]*track->fDz[ilayer])/(track->fSigmaZ[ilayer]*track->fSigmaZ[ilayer]) ; + // + for (cindex=0;cindex<100;cindex++){ + if (clusterindex[ilayer][cindex]<0) break; + if (clusterindex[ilayer][cindex]==tindex) break; + } + if (cindex>100) continue; + if (clusterweight[ilayer][cindex]>0.00001){ + sharefactor+= clusterweight[ilayer][cindex]; + cchi2/=clusterweight[ilayer][cindex]; + norm++; + } + newchi2+=cchi2; + } + newchi2/=(norm-track->fNSkipped-track->fNUsed); + track->fChi2MIP[4] = newchi2; + //if (norm>0) sharefactor/=norm; + //if (sharefactor<0.5) return 0; + // chi2*=1./(0.5+sharefactor); + if (newchi2GetNumberOfClusters()<4) { + return 0; + } + // + if ((minchi2/secchi2)<0.7){ + // + //increase the weight for clusters if the probability of other hypothesys is small + Double_t deltad = besttrack->GetD(GetX(),GetY()); + Double_t deltaz = besttrack->GetZat(GetX()) - GetZ(); + Double_t deltaprim = TMath::Sqrt(deltad*deltad+deltaz*deltaz); + // + for (Int_t icluster=0;iclusterGetNumberOfClusters();icluster++){ + Int_t index = besttrack->GetClusterIndex(icluster); + Int_t ilayer = (index & 0xf0000000) >> 28; + AliITSclusterV2 *c = (AliITSclusterV2*)GetCluster(index); + if (!c) continue; + if (c->GetNy()>3) continue; + if ( (ilayer*0.2+0.2)SetLabel(original->GetLabel()); + CookLabel(besttrack,0.0); + + // + return besttrack; +} +