//-------------------------------------------------------------------------
// Implementation of the ITS tracker class
-// It reads AliITSclusterV2 clusters and creates AliITStrackV2 tracks
+// It reads AliITSclusterV2 clusters and creates AliITStrackMI tracks
// and fills with them the ESD
-// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
-// Marian Ivanov, CERN, Marian.Ivanov@cern.ch
+// Origin: Marian Ivanov, CERN, Marian.Ivanov@cern.ch
// dEdx analysis by: Boris Batyunya, JINR, Boris.Batiounia@cern.ch
//
//-------------------------------------------------------------------------
#include "AliITSrecoV2.h"
#include <TTree.h>
#include "AliITSgeom.h"
-#include "AliTPCtrack.h"
#include "AliESD.h"
#include "AliITSclusterV2.h"
#include "AliITStrackerMI.h"
#include "TMatrixD.h"
+#include "TFile.h"
+#include "TTree.h"
#include "AliHelix.h"
-
-
+#include "AliESDV0MI.h"
+#include "AliLog.h"
+#include "TTreeStream.h"
ClassImp(AliITStrackerMI)
-ClassImp(AliITSRecV0Info)
//--------------------------------------------------------------------
//This is the AliITStrackerMI constructor
//--------------------------------------------------------------------
+ fCoeficients = 0;
+ fAfterV0 = kFALSE;
AliITSgeom *g=(AliITSgeom*)geom;
-
Float_t x,y,z;
Int_t i;
for (i=1; i<kMaxLayer+1; i++) {
for (Int_t i=0; i<kMaxLayer; i++) fLayersNotToSkip[i]=kLayersNotToSkip[i];
fLastLayerToTrackTo=kLastLayerToTrackTo;
+ for (Int_t i=0;i<100000;i++){
+ fBestTrackIndex[i]=0;
+ }
+ //
+ fDebugStreamer = new TTreeSRedirector("ITSdebug.root");
+
+}
+AliITStrackerMI::~AliITStrackerMI()
+{
+ //
+ //destructor
+ //
+ if (fCoeficients) delete []fCoeficients;
+ if (fDebugStreamer) {
+ //fDebugStreamer->Close();
+ delete fDebugStreamer;
+ }
}
void AliITStrackerMI::SetLayersNotToSkip(Int_t *l) {
for (Int_t i=0; i<kMaxLayer; i++) fgLayers[i].ResetClusters();
}
-static Int_t CorrectForDeadZoneMaterial(AliITStrackV2 *t) {
+static Int_t CorrectForDeadZoneMaterial(AliITStrackMI *t) {
//--------------------------------------------------------------------
// Correction for the material between the TPC and the ITS
// (should it belong to the TPC code ?)
// The clusters must be already loaded !
//--------------------------------------------------------------------
TObjArray itsTracks(15000);
-
+ fOriginal.Clear();
+ fEsd = event; // store pointer to the esd
{/* Read ESD tracks */
Int_t nentr=event->GetNumberOfTracks();
Info("Clusters2Tracks", "Number of ESD tracks: %d\n", nentr);
if ((esd->GetStatus()&AliESDtrack::kTPCin)==0) continue;
if (esd->GetStatus()&AliESDtrack::kTPCout) continue;
if (esd->GetStatus()&AliESDtrack::kITSin) continue;
-
- AliITStrackV2 *t=0;
+ if (esd->GetKinkIndex(0)>0) continue; //kink daughter
+ AliITStrackMI *t=0;
try {
- t=new AliITStrackV2(*esd);
+ t=new AliITStrackMI(*esd);
} catch (const Char_t *msg) {
- Warning("Clusters2Tracks",msg);
+ //Warning("Clusters2Tracks",msg);
delete t;
continue;
}
// write expected q
t->fExpQ = TMath::Max(0.8*t->fESDtrack->GetTPCsignal(),30.);
- if (TMath::Abs(t->fD[0])>10) {
- delete t;
- continue;
- }
-
- if (TMath::Abs(vdist)>20) {
- delete t;
- continue;
+ if (esd->GetV0Index(0)>0 && t->fD[0]<30){
+ //track - can be V0 according to TPC
}
- if (TMath::Abs(1/t->Get1Pt())<0.120) {
- delete t;
- continue;
- }
-
- if (CorrectForDeadZoneMaterial(t)!=0) {
- Warning("Clusters2Tracks",
- "failed to correct for the material in the dead zone !\n");
- delete t;
- continue;
+ else{
+ if (TMath::Abs(t->fD[0])>10) {
+ delete t;
+ continue;
+ }
+
+ if (TMath::Abs(vdist)>20) {
+ delete t;
+ continue;
+ }
+ if (TMath::Abs(1/t->Get1Pt())<0.120) {
+ delete t;
+ continue;
+ }
+
+ 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);
+ fOriginal.AddLast(t);
}
} /* End Read ESD tracks */
itsTracks.Sort();
+ fOriginal.Sort();
Int_t nentr=itsTracks.GetEntriesFast();
fTrackHypothesys.Expand(nentr);
+ fBestHypothesys.Expand(nentr);
MakeCoeficients(nentr);
Int_t ntrk=0;
for (fPass=0; fPass<2; fPass++) {
for (Int_t i=0; i<nentr; i++) {
// cerr<<fPass<<" "<<i<<'\r';
fCurrentEsdTrack = i;
- AliITStrackV2 *t=(AliITStrackV2*)itsTracks.UncheckedAt(i);
+ AliITStrackMI *t=(AliITStrackMI*)itsTracks.UncheckedAt(i);
if (t==0) continue; //this track has been already tracked
if (t->fReconstructed&&(t->fNUsed<1.5)) continue; //this track was already "succesfully" reconstructed
if ( (TMath::Abs(t->GetD(GetX(),GetY())) >3.) && fConstraint[fPass]) continue;
fI = 6;
ResetTrackToFollow(*t);
ResetBestTrack();
-
- FollowProlongationTree(t,i);
-
+ FollowProlongationTree(t,i,fConstraint[fPass]);
SortTrackHypothesys(fCurrentEsdTrack,20,0); //MI change
//
- AliITStrackV2 * besttrack = GetBestHypothesys(fCurrentEsdTrack,t,15);
+ AliITStrackMI * besttrack = GetBestHypothesys(fCurrentEsdTrack,t,15);
if (!besttrack) continue;
besttrack->SetLabel(tpcLabel);
// besttrack->CookdEdx();
CookdEdx(besttrack);
besttrack->fFakeRatio=1.;
CookLabel(besttrack,0.); //For comparison only
- // besttrack->UpdateESDtrack(AliESDtrack::kITSin);
- //
-
UpdateESDtrack(besttrack,AliESDtrack::kITSin);
-
- if ( besttrack->GetNumberOfClusters()<5 && fConstraint[fPass]) {
+
+ /*
+ if ( besttrack->GetNumberOfClusters()<6 && fConstraint[fPass]) {
continue;
}
if (besttrack->fChi2MIP[0]+besttrack->fNUsed>3.5) continue;
if ( (TMath::Abs(besttrack->fD[0]*besttrack->fD[0]+besttrack->fD[1]*besttrack->fD[1])>0.1) && fConstraint[fPass]) continue;
//delete itsTracks.RemoveAt(i);
+ */
+ if (fConstraint[fPass]&&(!besttrack->IsGoldPrimary())) continue; //to be tracked also without vertex constrain
+
+
t->fReconstructed = kTRUE;
ntrk++;
}
}
//GetBestHypothesysMIP(itsTracks);
- //FindV0(event);
-
+ UpdateTPCV0(event);
+ FindV02(event);
+ fAfterV0 = kTRUE;
+ //GetBestHypothesysMIP(itsTracks);
+ //
itsTracks.Delete();
//
Int_t entries = fTrackHypothesys.GetEntriesFast();
}
fTrackHypothesys.Delete();
+ fBestHypothesys.Delete();
+ fOriginal.Clear();
+ delete []fCoeficients;
+ fCoeficients=0;
Info("Clusters2Tracks","Number of prolonged tracks: %d\n",ntrk);
-
+
return 0;
}
-
-Int_t AliITStrackerMI::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 */
- 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();
-
- Info("Clusters2Tracks","Number of TPC tracks: %d\n",nentr);
-
- for (Int_t i=0; i<nentr; i++) {
- tpcTree->GetEvent(i);
- AliITStrackV2 *t=0;
- try {
- t=new AliITStrackV2(*itrack);
- } catch (const Char_t *msg) {
- Warning("Clusters2Tracks",msg);
- delete t;
- 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 itrack;
- }
- itsTracks.Sort();
- nentr=itsTracks.GetEntriesFast();
-
-
- AliITStrackV2 *otrack=&fBestTrack;
- 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;
- for (Int_t i=0; i<nentr; i++) {
- AliITStrackV2 *t=(AliITStrackV2*)itsTracks.UncheckedAt(i);
- if (t==0) continue; //this track has been already tracked
- Int_t tpcLabel=t->GetLabel(); //save the TPC track label
-
- ResetTrackToFollow(*t);
- ResetBestTrack();
- /*
- for (FollowProlongation(); fI<kMaxLayer; fI++) {
- while (TakeNextProlongation()) FollowProlongation();
- }
- */
- FollowProlongationTree(t,i);
- if (fBestTrack.GetNumberOfClusters() == 0) continue;
-
- if (fConstraint[fPass]) {
- ResetTrackToFollow(*t);
- if (!RefitAt(3.7, &fTrackToFollow, &fBestTrack)) continue;
- ResetBestTrack();
- }
-
- fBestTrack.SetLabel(tpcLabel);
- //fBestTrack.CookdEdx();
- CookdEdx(&fBestTrack);
-
- CookLabel(&fBestTrack,0.); //For comparison only
- itsTree->Fill();
- //UseClusters(&fBestTrack);
- delete itsTracks.RemoveAt(i);
- }
- }
-
- nentr=(Int_t)itsTree->GetEntries();
- Info("Clusters2Tracks","Number of prolonged tracks: %d\n",nentr);
-
- itsTracks.Delete();
-
- return 0;
-}
-
Int_t AliITStrackerMI::PropagateBack(AliESD *event) {
//--------------------------------------------------------------------
// This functions propagates reconstructed ITS tracks back
if ((esd->GetStatus()&AliESDtrack::kITSin)==0) continue;
if (esd->GetStatus()&AliESDtrack::kITSout) continue;
- AliITStrackV2 *t=0;
+ AliITStrackMI *t=0;
try {
- t=new AliITStrackV2(*esd);
+ t=new AliITStrackMI(*esd);
} catch (const Char_t *msg) {
- Warning("PropagateBack",msg);
+ //Warning("PropagateBack",msg);
delete t;
continue;
}
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");
+ //Warning("PropagateBack",
+ // "failed to correct for the material in the dead zone !\n");
delete t;
continue;
}
if (esd->GetStatus()&AliESDtrack::kTPCout)
if ((esd->GetStatus()&AliESDtrack::kTPCrefit)==0) continue;
- AliITStrackV2 *t=0;
+ AliITStrackMI *t=0;
try {
- t=new AliITStrackV2(*esd);
+ t=new AliITStrackMI(*esd);
} catch (const Char_t *msg) {
- Warning("RefitInward",msg);
+ //Warning("RefitInward",msg);
delete t;
continue;
}
t->fExpQ = TMath::Max(0.8*t->fESDtrack->GetTPCsignal(),30.);
if (CorrectForDeadZoneMaterial(t)!=0) {
- Warning("RefitInward",
- "failed to correct for the material in the dead zone !\n");
+ //Warning("RefitInward",
+ // "failed to correct for the material in the dead zone !\n");
delete t;
continue;
}
}
-void AliITStrackerMI::FollowProlongationTree(AliITStrackV2 * otrack, Int_t esdindex)
+void AliITStrackerMI::FollowProlongationTree(AliITStrackMI * otrack, Int_t esdindex, Bool_t constrain)
{
//--------------------------------------------------------------------
// Follow prolongation tree
//--------------------------------------------------------------------
+ //
+ AliESDtrack * esd = otrack->fESDtrack;
+ if (esd->GetV0Index(0)>0){
+ //
+ // TEMPORARY SOLLUTION: map V0 indexes to point to proper track
+ // mapping of esd track is different as its track in Containers
+ // Need something more stable
+ // Indexes are set back againg to the ESD track indexes in UpdateTPCV0
+ for (Int_t i=0;i<3;i++){
+ Int_t index = esd->GetV0Index(i);
+ if (index==0) break;
+ AliESDV0MI * vertex = fEsd->GetV0MI(index);
+ if (vertex->GetStatus()<0) continue; // rejected V0
+ //
+ if (esd->GetSign()>0) {
+ vertex->SetIndex(0,esdindex);
+ }
+ else{
+ vertex->SetIndex(1,esdindex);
+ }
+ }
+ }
+ TObjArray *bestarray = (TObjArray*)fBestHypothesys.At(esdindex);
+ if (!bestarray){
+ bestarray = new TObjArray(5);
+ fBestHypothesys.AddAt(bestarray,esdindex);
+ }
+ //
//setup tree of the prolongations
//
- static AliITStrackV2 tracks[7][100];
- AliITStrackV2 *currenttrack;
- static AliITStrackV2 currenttrack1;
- static AliITStrackV2 currenttrack2;
- static AliITStrackV2 backuptrack;
+ static AliITStrackMI tracks[7][100];
+ AliITStrackMI *currenttrack;
+ static AliITStrackMI currenttrack1;
+ static AliITStrackMI currenttrack2;
+ static AliITStrackMI backuptrack;
Int_t ntracks[7];
Int_t nindexes[7][100];
Float_t normalizedchi2[100];
for (Int_t ilayer=0;ilayer<6;ilayer++) ntracks[ilayer]=0;
otrack->fNSkipped=0;
- new (&(tracks[6][0])) AliITStrackV2(*otrack);
+ new (&(tracks[6][0])) AliITStrackMI(*otrack);
ntracks[6]=1;
- nindexes[6][0]=0;
+ for (Int_t i=0;i<7;i++) nindexes[i][0]=0;
//
//
// follow prolongations
if (itrack>1&&tracks[ilayer+1][nindexes[ilayer+1][itrack]].fNUsed>2. && nused>3) continue;
}
- new(¤ttrack1) AliITStrackV2(tracks[ilayer+1][nindexes[ilayer+1][itrack]]);
+ new(¤ttrack1) AliITStrackMI(tracks[ilayer+1][nindexes[ilayer+1][itrack]]);
if (ilayer==3 || ilayer==1) {
Double_t rs=0.5*(fgLayers[ilayer+1].GetR() + r);
Double_t d=0.0034, x0=38.6;
//propagate to the intersection
const AliITSdetector &det=layer.GetDetector(idet);
phi=det.GetPhi();
- new(¤ttrack2) AliITStrackV2(currenttrack1);
+ new(¤ttrack2) AliITStrackMI(currenttrack1);
if (!currenttrack1.Propagate(phi,det.GetR())) {
continue;
}
//
Double_t msz=1./((currenttrack1.GetSigmaZ2() + 16.*kSigmaZ2[ilayer]));
Double_t msy=1./((currenttrack1.GetSigmaY2() + 16.*kSigmaY2[ilayer]));
- if (fConstraint[fPass]){
+ if (constrain){
msy/=60; msz/=60.;
}
else{
Float_t pz = (z - c->GetZ()) , py=(y - c->GetY());
if (pz*pz*msz+py*py*msy>1.) continue;
//
- new (&backuptrack) AliITStrackV2(currenttrack2);
+ new (&backuptrack) AliITStrackMI(currenttrack2);
change = kTRUE;
currenttrack =¤ttrack2;
if (!currenttrack->Propagate(det.GetPhi(),det.GetR())) {
- new (currenttrack) AliITStrackV2(backuptrack);
+ new (currenttrack) AliITStrackMI(backuptrack);
change = kFALSE;
continue;
}
if (chi2<kMaxChi2s[ilayer]){
if (c->GetQ()==0) deadzone=1; // take dead zone only once
if (ntracks[ilayer]>=100) continue;
- AliITStrackV2 * updatetrack = new (&tracks[ilayer][ntracks[ilayer]]) AliITStrackV2(*currenttrack);
+ AliITStrackMI * updatetrack = new (&tracks[ilayer][ntracks[ilayer]]) AliITStrackMI(*currenttrack);
updatetrack->fClIndex[ilayer]=0;
if (change){
- new (¤ttrack2) AliITStrackV2(backuptrack);
+ new (¤ttrack2) AliITStrackMI(backuptrack);
}
if (c->GetQ()!=0){
if (!UpdateMI(updatetrack,c,chi2,(ilayer<<28)+ci)) continue;
Double_t x0;
Double_t d=layer.GetThickness(updatetrack->GetY(),updatetrack->GetZ(),x0);
updatetrack->CorrectForMaterial(d,x0);
- if (fConstraint[fPass]) {
- updatetrack->fConstrain = fConstraint[fPass];
+ if (constrain) {
+ updatetrack->fConstrain = constrain;
fI = ilayer;
Double_t d=GetEffectiveThickness(0,0); //Think of this !!!!
Double_t xyz[]={GetX(),GetY(),GetZ()};
ntracks[ilayer]++;
} // create new hypothesy
} // loop over possible cluster prolongation
- // if (fConstraint[fPass]&&itrack<2&¤ttrack1.fNSkipped==0 && deadzone==0){
- if (itrack<2&¤ttrack1.fNSkipped==0 && deadzone==0&&ntracks[ilayer]<100){
- AliITStrackV2* vtrack = new (&tracks[ilayer][ntracks[ilayer]]) AliITStrackV2(currenttrack1);
+ // if (constrain&&itrack<2&¤ttrack1.fNSkipped==0 && deadzone==0){
+ if (constrain&&itrack<2&¤ttrack1.fNSkipped==0 && deadzone==0&&ntracks[ilayer]<100){
+ AliITStrackMI* vtrack = new (&tracks[ilayer][ntracks[ilayer]]) AliITStrackMI(currenttrack1);
vtrack->fClIndex[ilayer]=0;
fI = ilayer;
Double_t d=GetEffectiveThickness(0,0); //Think of this !!!!
vtrack->fNSkipped++;
ntracks[ilayer]++;
}
+
+ if (constrain&&itrack<1&&TMath::Abs(currenttrack1.fP3)>1.1){ //big theta -- for low mult. runs
+ AliITStrackMI* vtrack = new (&tracks[ilayer][ntracks[ilayer]]) AliITStrackMI(currenttrack1);
+ vtrack->fClIndex[ilayer]=0;
+ fI = ilayer;
+ Double_t d=GetEffectiveThickness(0,0); //Think of this !!!!
+ Double_t xyz[]={GetX(),GetY(),GetZ()};
+ Double_t ers[]={GetSigmaX(),GetSigmaY(),GetSigmaZ()};
+ vtrack->Improve(d,xyz,ers);
+ vtrack->fNDeadZone++;
+ ntracks[ilayer]++;
+ }
+
} //loop over track candidates
//
if ( normalizedchi2[itrack]<3+0.5*ilayer) golds++;
if (ilayer>4) accepted++;
else{
- if ( fConstraint[fPass] && normalizedchi2[itrack]<kMaxNormChi2C[ilayer]+1) accepted++;
- if (!fConstraint[fPass] && normalizedchi2[itrack]<kMaxNormChi2NonC[ilayer]+1) accepted++;
+ if ( constrain && normalizedchi2[itrack]<kMaxNormChi2C[ilayer]+1) accepted++;
+ if (!constrain && normalizedchi2[itrack]<kMaxNormChi2NonC[ilayer]+1) accepted++;
}
}
TMath::Sort(ntracks[ilayer],normalizedchi2,nindexes[ilayer],kFALSE);
if (ntracks[ilayer]>90) ntracks[ilayer]=90;
} //loop over layers
//printf("%d\t%d\t%d\t%d\t%d\t%d\n",ntracks[0],ntracks[1],ntracks[2],ntracks[3],ntracks[4],ntracks[5]);
+ Int_t max = constrain? 20: 5;
- for (Int_t i=0;i<TMath::Min(20,ntracks[0]);i++) {
- AliITStrackV2 & track= tracks[0][nindexes[0][i]];
- if (!fConstraint[fPass]&&track.fNormChi2[0]>7.)continue;
- AddTrackHypothesys(new AliITStrackV2(track), esdindex);
+ for (Int_t i=0;i<TMath::Min(max,ntracks[0]);i++) {
+ AliITStrackMI & track= tracks[0][nindexes[0][i]];
+ if (track.GetNumberOfClusters()<2) continue;
+ if (!constrain&&track.fNormChi2[0]>7.)continue;
+ AddTrackHypothesys(new AliITStrackMI(track), esdindex);
}
- for (Int_t i=0;i<TMath::Min(4,ntracks[1]);i++) {
- AliITStrackV2 & track= tracks[1][nindexes[1][i]];
+ for (Int_t i=0;i<TMath::Min(2,ntracks[1]);i++) {
+ AliITStrackMI & track= tracks[1][nindexes[1][i]];
if (track.GetNumberOfClusters()<4) continue;
- if (!fConstraint[fPass]&&track.fNormChi2[1]>7.)continue;
- if (fConstraint[fPass]) track.fNSkipped+=1;
- if (!fConstraint[fPass]) {
+ if (!constrain&&track.fNormChi2[1]>7.)continue;
+ if (constrain) track.fNSkipped+=1;
+ if (!constrain) {
track.fD[0] = track.GetD(GetX(),GetY());
track.fNSkipped+=4./(4.+8.*TMath::Abs(track.fD[0]));
if (track.fN+track.fNDeadZone+track.fNSkipped>6) {
track.fNSkipped = 6-track.fN+track.fNDeadZone;
}
}
- AddTrackHypothesys(new AliITStrackV2(track), esdindex);
+ AddTrackHypothesys(new AliITStrackMI(track), esdindex);
}
//}
- if (!fConstraint[fPass]){
- for (Int_t i=0;i<TMath::Min(3,ntracks[2]);i++) {
- AliITStrackV2 & track= tracks[2][nindexes[2][i]];
- if (track.GetNumberOfClusters()<4) continue;
- if (!fConstraint[fPass]&&track.fNormChi2[2]>7.)continue;
- if (fConstraint[fPass]) track.fNSkipped+=2;
- if (!fConstraint[fPass]){
+ if (!constrain){
+ for (Int_t i=0;i<TMath::Min(2,ntracks[2]);i++) {
+ AliITStrackMI & track= tracks[2][nindexes[2][i]];
+ if (track.GetNumberOfClusters()<3) continue;
+ if (!constrain&&track.fNormChi2[2]>7.)continue;
+ if (constrain) track.fNSkipped+=2;
+ if (!constrain){
track.fD[0] = track.GetD(GetX(),GetY());
track.fNSkipped+= 7./(7.+8.*TMath::Abs(track.fD[0]));
if (track.fN+track.fNDeadZone+track.fNSkipped>6) {
track.fNSkipped = 6-track.fN+track.fNDeadZone;
}
}
- AddTrackHypothesys(new AliITStrackV2(track), esdindex);
+ AddTrackHypothesys(new AliITStrackMI(track), esdindex);
}
}
+
+ if (!constrain){
+ //
+ // register best tracks - important for V0 finder
+ //
+ for (Int_t ilayer=0;ilayer<5;ilayer++){
+ if (ntracks[ilayer]==0) continue;
+ AliITStrackMI & track= tracks[ilayer][nindexes[ilayer][0]];
+ if (track.GetNumberOfClusters()<1) continue;
+ CookLabel(&track,0);
+ bestarray->AddAt(new AliITStrackMI(track),ilayer);
+ }
+ }
+ //
+ // update TPC V0 information
+ //
+ if (otrack->fESDtrack->GetV0Index(0)>0){
+ Float_t fprimvertex[3]={GetX(),GetY(),GetZ()};
+ for (Int_t i=0;i<3;i++){
+ Int_t index = otrack->fESDtrack->GetV0Index(i);
+ if (index==0) break;
+ AliESDV0MI * vertex = fEsd->GetV0MI(index);
+ if (vertex->GetStatus()<0) continue; // rejected V0
+ //
+ if (otrack->fP4>0) {
+ vertex->SetIndex(0,esdindex);
+ }
+ else{
+ vertex->SetIndex(1,esdindex);
+ }
+ //find nearest layer with track info
+ Int_t nearestold = GetNearestLayer(vertex->GetXrp());
+ Int_t nearest = nearestold;
+ for (Int_t ilayer =nearest;ilayer<8;ilayer++){
+ if (ntracks[nearest]==0){
+ nearest = ilayer;
+ }
+ }
+ //
+ AliITStrackMI & track= tracks[nearest][nindexes[nearest][0]];
+ if (nearestold<5&&nearest<5){
+ Bool_t accept = track.fNormChi2[nearest]<10;
+ if (accept){
+ if (track.fP4>0) {
+ vertex->SetP(track);
+ vertex->Update(fprimvertex);
+ // vertex->SetIndex(0,track.fESDtrack->GetID());
+ if (track.GetNumberOfClusters()>2) AddTrackHypothesys(new AliITStrackMI(track), esdindex);
+ }else{
+ vertex->SetM(track);
+ vertex->Update(fprimvertex);
+ //vertex->SetIndex(1,track.fESDtrack->GetID());
+ if (track.GetNumberOfClusters()>2) AddTrackHypothesys(new AliITStrackMI(track), esdindex);
+ }
+ vertex->SetStatus(vertex->GetStatus()+1);
+ }else{
+ // vertex->SetStatus(-2); // reject V0 - not enough clusters
+ }
+ }
+ // if (nearestold>3){
+// Int_t indexlayer = (ntracks[0]>0)? 0:1;
+// if (ntracks[indexlayer]>0){
+// AliITStrackMI & track= tracks[indexlayer][nindexes[indexlayer][0]];
+// if (track.GetNumberOfClusters()>4&&track.fNormChi2[indexlayer]<4){
+// vertex->SetStatus(-1); // reject V0 - clusters before
+// }
+// }
+// }
+ }
+ }
}
if (n>1) fRoad=2*fR*TMath::Sqrt(3.14/n);
}
+
Int_t AliITStrackerMI::AliITSlayer::InsertCluster(AliITSclusterV2 *c) {
//--------------------------------------------------------------------
//This function adds a cluster to this layer
return 1;
}
fCurrentSlice=-1;
- if (fN==0) {fClusters[fN++]=c; return 0;}
- Int_t i=FindClusterIndex(c->GetZ());
- memmove(fClusters+i+1 ,fClusters+i,(fN-i)*sizeof(AliITSclusterV2*));
- memmove(fY+i+1 ,fY+i,(fN-i)*sizeof(Float_t));
- memmove(fZ+i+1 ,fZ+i,(fN-i)*sizeof(Float_t));
- fClusters[i]=c; fN++;
- //
+ fClusters[fN]=c;
+ fN++;
AliITSdetector &det=GetDetector(c->GetDetectorIndex());
- Double_t y=fR*det.GetPhi() + c->GetY();
- if (y>2.*fR*TMath::Pi()) y -= 2.*fR*TMath::Pi();
- fY[i] = y;
- fZ[i] = c->GetZ();
if (c->GetY()<det.GetYmin()) det.SetYmin(c->GetY());
if (c->GetY()>det.GetYmax()) det.SetYmax(c->GetY());
if (c->GetZ()<det.GetZmin()) det.SetZmin(c->GetZ());
//
//sort clusters
//
+ AliITSclusterV2 **clusters = new AliITSclusterV2*[fN];
+ Float_t *z = new Float_t[fN];
+ Int_t * index = new Int_t[fN];
+ //
+ for (Int_t i=0;i<fN;i++){
+ z[i] = fClusters[i]->GetZ();
+ }
+ TMath::Sort(fN,z,index,kFALSE);
+ for (Int_t i=0;i<fN;i++){
+ clusters[i] = fClusters[index[i]];
+ }
+ //
+ for (Int_t i=0;i<fN;i++){
+ fClusters[i] = clusters[i];
+ fZ[i] = fClusters[i]->GetZ();
+ AliITSdetector &det=GetDetector(fClusters[i]->GetDetectorIndex());
+ Double_t y=fR*det.GetPhi() + fClusters[i]->GetY();
+ if (y>2.*fR*TMath::Pi()) y -= 2.*fR*TMath::Pi();
+ fY[i] = y;
+ }
+ delete[] index;
+ delete[] z;
+ delete[] clusters;
+ //
+
fYB[0]=10000000;
fYB[1]=-10000000;
for (Int_t i=0;i<fN;i++){
for (Int_t i=0;i<21;i++) {fBy20[i][0] = fYB[0]+(i-0.75)*fDy20; fBy20[i][1] = fYB[0]+(i+0.75)*fDy20;}
//
//
- for (Int_t i=0;i<fN;i++) for (Int_t irot=-1;irot<=1;irot++){
- Float_t curY = fY[i]+irot*2.*TMath::Pi()*fR;
- if (curY<fYB[0]-fDy5) continue;
- if (curY>fYB[1]+fDy5) continue;
- //
- // slice 10
- if (TMath::Abs(fYB[1]-fYB[0])<=0) continue;
- Float_t fslice = TMath::Nint(10.*(curY-fYB[0])/(fYB[1]-fYB[0]));
- Float_t ymiddle = fYB[0]+fslice*fDy10;
- for (Int_t di =-1;di<=1;di++){
- if (TMath::Abs(curY-(ymiddle+(float)di*fDy10))<0.75*fDy10){
- //
- Int_t slice = int(fslice+21.0001)-21+di;
- if (slice<0) continue;
- if (slice>10) continue;
- fClusters10[slice][fN10[slice]] = fClusters[i];
- fY10[slice][fN10[slice]] = curY;
- fZ10[slice][fN10[slice]] = fZ[i];
- fClusterIndex10[slice][fN10[slice]]=i;
- fN10[slice]++;
+ for (Int_t i=0;i<fN;i++)
+ for (Int_t irot=-1;irot<=1;irot++){
+ Float_t curY = fY[i]+irot*TMath::TwoPi()*fR;
+ // slice 5
+ for (Int_t slice=0; slice<6;slice++){
+ if (fBy5[slice][0]<curY && curY<fBy5[slice][1]&&fN5[slice]<kMaxClusterPerLayer5){
+ fClusters5[slice][fN5[slice]] = fClusters[i];
+ fY5[slice][fN5[slice]] = curY;
+ fZ5[slice][fN5[slice]] = fZ[i];
+ fClusterIndex5[slice][fN5[slice]]=i;
+ fN5[slice]++;
+ }
}
- }
- //
- // slice 5
- fslice = TMath::Nint(5.*(curY-fYB[0])/(fYB[1]-fYB[0]));
- ymiddle = fYB[0]+fslice*fDy5;
- for (Int_t di =-1;di<=1;di++){
- if (TMath::Abs(curY-(ymiddle+(float)di*fDy5))<0.75*fDy5){
- //
- Int_t slice = int(fslice+21.0001)-21+di;
- if (slice<0) continue;
- if (slice>5) continue;
- fClusters5[slice][fN5[slice]] = fClusters[i];
- fY5[slice][fN5[slice]] = curY;
- fZ5[slice][fN5[slice]] = fZ[i];
- fClusterIndex5[slice][fN5[slice]]=i;
- fN5[slice]++;
+ // slice 10
+ for (Int_t slice=0; slice<11;slice++){
+ if (fBy10[slice][0]<curY && curY<fBy10[slice][1]&&fN10[slice]<kMaxClusterPerLayer10){
+ fClusters10[slice][fN10[slice]] = fClusters[i];
+ fY10[slice][fN10[slice]] = curY;
+ fZ10[slice][fN10[slice]] = fZ[i];
+ fClusterIndex10[slice][fN10[slice]]=i;
+ fN10[slice]++;
+ }
}
+ // slice 20
+ for (Int_t slice=0; slice<21;slice++){
+ if (fBy20[slice][0]<curY && curY<fBy20[slice][1]&&fN20[slice]<kMaxClusterPerLayer20){
+ fClusters20[slice][fN20[slice]] = fClusters[i];
+ fY20[slice][fN20[slice]] = curY;
+ fZ20[slice][fN20[slice]] = fZ[i];
+ fClusterIndex20[slice][fN20[slice]]=i;
+ fN20[slice]++;
+ }
+ }
}
- //
- // slice 20
- fslice = TMath::Nint(20.*(curY-fYB[0])/(fYB[1]-fYB[0]));
- ymiddle = fYB[0]+fslice*fDy20;
- for (Int_t di =-1;di<=1;di++){
- if (TMath::Abs(curY-(ymiddle+(float)di*fDy20))<0.75*fDy20){
- //
- Int_t slice = int(fslice+21.0001)-21+di;
- if (slice<0) continue;
- if (slice>20) continue;
- fClusters20[slice][fN20[slice]] = fClusters[i];
- fY20[slice][fN20[slice]] = curY;
- fZ20[slice][fN20[slice]] = fZ[i];
- fClusterIndex20[slice][fN20[slice]]=i;
- fN20[slice]++;
- }
+
+ //
+ // consistency check
+ //
+ for (Int_t i=0;i<fN-1;i++){
+ if (fZ[i]>fZ[i+1]){
+ printf("Bugg\n");
+ }
+ }
+ //
+ for (Int_t slice=0;slice<21;slice++)
+ for (Int_t i=0;i<fN20[slice]-1;i++){
+ if (fZ20[slice][i]>fZ20[slice][i+1]){
+ printf("Bugg\n");
}
}
-}
+}
Int_t AliITStrackerMI::AliITSlayer::FindClusterIndex(Float_t z) const {
}
return m;
}
-/*
-void AliITStrackerMI::AliITSlayer::
-SelectClusters(Double_t zmin,Double_t zmax,Double_t ymin, Double_t ymax) {
- //--------------------------------------------------------------------
- // This function sets the "window"
- //--------------------------------------------------------------------
- fI=FindClusterIndex(zmin); fZmax=zmax;
- fImax = TMath::Min(FindClusterIndex(zmax)+1,fN);
- Double_t circle=2*TMath::Pi()*fR;
- if (ymax>circle) { ymax-=circle; ymin-=circle; }
- fYmin=ymin; fYmax=ymax;
- fSkip = 0;
-}
-*/
void AliITStrackerMI::AliITSlayer::
fAccepted =0;
}
-/*
-const AliITSclusterV2 *AliITStrackerMI::AliITSlayer::GetNextCluster(Int_t &ci){
+
+
+
+Int_t AliITStrackerMI::AliITSlayer::
+FindDetectorIndex(Double_t phi, Double_t z) const {
//--------------------------------------------------------------------
- // This function returns clusters within the "window"
+ //This function finds the detector crossed by the track
//--------------------------------------------------------------------
- const AliITSclusterV2 *cluster=0;
- for (Int_t i=fI; i<fN; i++) {
- const AliITSclusterV2 *c=fClusters[i];
- if (c->GetZ() > fZmax) break;
- // if (c->IsUsed()) continue;
- const AliITSdetector &det=GetDetector(c->GetDetectorIndex());
- Double_t y=fR*det.GetPhi() + c->GetY();
+ 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);
+ if (np>=fNladders) np-=fNladders;
+ if (np<0) np+=fNladders;
- if (y>2.*fR*TMath::Pi()) y -= 2*fR*TMath::Pi();
- if (y>1.*fR*TMath::Pi() && fYmax<y) y -= 2*fR*TMath::Pi();
+ Double_t dz=fZOffset-z;
+ Int_t nz=Int_t(dz*(fNdetectors-1)*0.5/fZOffset+0.5);
+ if (nz>=fNdetectors) return -1;
+ if (nz<0) return -1;
- if (y<fYmin) continue;
- if (y>fYmax) continue;
- cluster=c; ci=i;
- fI=i+1;
- break;
- }
- return cluster;
+ return np*fNdetectors + nz;
}
-*/
+
const AliITSclusterV2 *AliITStrackerMI::AliITSlayer::GetNextCluster(Int_t &ci){
//--------------------------------------------------------------------
for (Int_t i=fI; i<fImax; i++) {
Double_t y = fY[i];
if (fYmax<y) y -= rpi2;
+ if (fYmin>y) y += rpi2;
if (y<fYmin) continue;
if (y>fYmax) continue;
if (fClusters[i]->GetQ()==0&&fSkip==2) continue;
-Int_t AliITStrackerMI::AliITSlayer::
-FindDetectorIndex(Double_t phi, Double_t z) const {
- //--------------------------------------------------------------------
- //This function finds the detector crossed by the track
- //--------------------------------------------------------------------
- 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);
- if (np>=fNladders) np-=fNladders;
- if (np<0) np+=fNladders;
-
- Double_t dz=fZOffset-z;
- Int_t nz=Int_t(dz*(fNdetectors-1)*0.5/fZOffset+0.5);
- if (nz>=fNdetectors) return -1;
- if (nz<0) return -1;
-
- return np*fNdetectors + nz;
-}
-
Double_t AliITStrackerMI::AliITSlayer::GetThickness(Double_t y,Double_t z,Double_t &x0)
const {
//--------------------------------------------------------------------
}
Bool_t
-AliITStrackerMI::RefitAt(Double_t xx,AliITStrackV2 *t,const AliITStrackV2 *c) {
+AliITStrackerMI::RefitAt(Double_t xx,AliITStrackMI *t,const AliITStrackMI *c) {
//--------------------------------------------------------------------
// This function refits the track "t" at the position "x" using
// the clusters from "c"
}
-Double_t AliITStrackerMI::GetNormalizedChi2(AliITStrackV2 * track, Int_t mode)
+Double_t AliITStrackerMI::GetNormalizedChi2(AliITStrackMI * track, Int_t mode)
{
//
// calculate normalized chi2
}
-Double_t AliITStrackerMI::GetMatchingChi2(AliITStrackV2 * track1, AliITStrackV2 * track2)
+Double_t AliITStrackerMI::GetMatchingChi2(AliITStrackMI * track1, AliITStrackMI * track2)
{
//
// return matching chi2 between two tracks
- AliITStrackV2 track3(*track2);
+ AliITStrackMI track3(*track2);
track3.Propagate(track1->GetAlpha(),track1->GetX());
TMatrixD vec(5,1);
vec(0,0)=track1->fP0-track3.fP0;
}
-Double_t AliITStrackerMI::GetTruncatedChi2(AliITStrackV2 * track, Float_t fac)
+Double_t AliITStrackerMI::GetTruncatedChi2(AliITStrackMI * track, Float_t fac)
{
//
// calculate normalized chi2
}
-Double_t AliITStrackerMI::GetInterpolatedChi2(AliITStrackV2 * forwardtrack, AliITStrackV2 * backtrack)
+Double_t AliITStrackerMI::GetInterpolatedChi2(AliITStrackMI * forwardtrack, AliITStrackMI * backtrack)
{
//
// calculate normalized chi2
return fgLayers[l].GetWeight(c);
}
-void AliITStrackerMI::RegisterClusterTracks(AliITStrackV2* track,Int_t id)
+void AliITStrackerMI::RegisterClusterTracks(AliITStrackMI* track,Int_t id)
{
//---------------------------------------------
// register track to the list
+ //
+ if (track->fESDtrack->GetKinkIndex(0)!=0) return; //don't register kink tracks
+ //
+ //
for (Int_t icluster=0;icluster<track->GetNumberOfClusters();icluster++){
Int_t index = track->GetClusterIndex(icluster);
Int_t l=(index & 0xf0000000) >> 28;
}
}
}
-void AliITStrackerMI::UnRegisterClusterTracks(AliITStrackV2* track, Int_t id)
+void AliITStrackerMI::UnRegisterClusterTracks(AliITStrackMI* track, Int_t id)
{
//---------------------------------------------
// unregister track from the list
}
}
}
-Float_t AliITStrackerMI::GetNumberOfSharedClusters(AliITStrackV2* track,Int_t id, Int_t list[6], AliITSclusterV2 *clist[6])
+Float_t AliITStrackerMI::GetNumberOfSharedClusters(AliITStrackMI* track,Int_t id, Int_t list[6], AliITSclusterV2 *clist[6])
{
//-------------------------------------------------------------
//get number of shared clusters
return shared;
}
-Int_t AliITStrackerMI::GetOverlapTrack(AliITStrackV2 *track, Int_t trackID, Int_t &shared, Int_t clusterlist[6],Int_t overlist[6])
+Int_t AliITStrackerMI::GetOverlapTrack(AliITStrackMI *track, Int_t trackID, Int_t &shared, Int_t clusterlist[6],Int_t overlist[6])
{
//
// find first shared track
}
-AliITStrackV2 * AliITStrackerMI::GetBest2Tracks(Int_t trackID1, Int_t trackID2, Float_t th0, Float_t th1){
+AliITStrackMI * AliITStrackerMI::GetBest2Tracks(Int_t trackID1, Int_t trackID2, Float_t th0, Float_t th1){
//
// try to find track hypothesys without conflicts
// with minimal chi2;
TClonesArray *arr1 = (TClonesArray*)fTrackHypothesys.At(trackID1);
Int_t entries1 = arr1->GetEntriesFast();
TClonesArray *arr2 = (TClonesArray*)fTrackHypothesys.At(trackID2);
- if (!arr2) return (AliITStrackV2*) arr1->UncheckedAt(0);
+ if (!arr2) return (AliITStrackMI*) arr1->UncheckedAt(0);
Int_t entries2 = arr2->GetEntriesFast();
- if (entries2<=0) return (AliITStrackV2*) arr1->UncheckedAt(0);
+ if (entries2<=0) return (AliITStrackMI*) arr1->UncheckedAt(0);
//
- AliITStrackV2 * track10=(AliITStrackV2*) arr1->UncheckedAt(0);
- AliITStrackV2 * track20=(AliITStrackV2*) arr2->UncheckedAt(0);
+ AliITStrackMI * track10=(AliITStrackMI*) arr1->UncheckedAt(0);
+ AliITStrackMI * track20=(AliITStrackMI*) arr2->UncheckedAt(0);
if (TMath::Abs(1./track10->Get1Pt())>0.5+TMath::Abs(1/track20->Get1Pt())) return track10;
for (Int_t itrack=0;itrack<entries1;itrack++){
- AliITStrackV2 * track=(AliITStrackV2*) arr1->UncheckedAt(itrack);
+ AliITStrackMI * track=(AliITStrackMI*) arr1->UncheckedAt(itrack);
UnRegisterClusterTracks(track,trackID1);
}
//
for (Int_t itrack=0;itrack<entries2;itrack++){
- AliITStrackV2 * track=(AliITStrackV2*) arr2->UncheckedAt(itrack);
+ AliITStrackMI * track=(AliITStrackMI*) arr2->UncheckedAt(itrack);
UnRegisterClusterTracks(track,trackID2);
}
Int_t index1=0;
Float_t * ny2 = GetNy(trackID2), * nz2 = GetNz(trackID2);
for (Int_t itrack1=0;itrack1<entries1;itrack1++){
- AliITStrackV2 * track1=(AliITStrackV2*) arr1->UncheckedAt(itrack1);
+ AliITStrackMI * track1=(AliITStrackMI*) arr1->UncheckedAt(itrack1);
//if (track1->fFakeRatio>0) continue;
RegisterClusterTracks(track1,trackID1);
for (Int_t itrack2=0;itrack2<entries2;itrack2++){
- AliITStrackV2 * track2=(AliITStrackV2*) arr2->UncheckedAt(itrack2);
+ AliITStrackMI * track2=(AliITStrackMI*) arr2->UncheckedAt(itrack2);
// Float_t current = w1*track1->fChi2MIP[0]+w2*track2->fChi2MIP[0];
//if (track2->fFakeRatio>0) continue;
// if (maxconflicts<4 && maxchi2<th0){
if (maxchi2<th0*2.){
Float_t orig = track10->fFakeRatio*track10->GetNumberOfClusters();
- AliITStrackV2* track1=(AliITStrackV2*) arr1->UncheckedAt(index1);
+ AliITStrackMI* track1=(AliITStrackMI*) arr1->UncheckedAt(index1);
track1->fChi2MIP[5] = maxconflicts;
track1->fChi2MIP[6] = maxchi2;
track1->fChi2MIP[7] = 0.01+orig-(track1->fFakeRatio*track1->GetNumberOfClusters());
}
for (Int_t itrack=0;itrack<entries1;itrack++){
- AliITStrackV2 * track=(AliITStrackV2*) arr1->UncheckedAt(itrack);
+ AliITStrackMI * track=(AliITStrackMI*) arr1->UncheckedAt(itrack);
UnRegisterClusterTracks(track,trackID1);
}
//
for (Int_t itrack=0;itrack<entries2;itrack++){
- AliITStrackV2 * track=(AliITStrackV2*) arr2->UncheckedAt(itrack);
+ AliITStrackMI * track=(AliITStrackMI*) arr2->UncheckedAt(itrack);
UnRegisterClusterTracks(track,trackID2);
}
}
-void AliITStrackerMI::AddTrackHypothesys(AliITStrackV2 * track, Int_t esdindex)
+void AliITStrackerMI::AddTrackHypothesys(AliITStrackMI * track, Int_t esdindex)
{
//------------------------------------------------------------------
// add track to the list of hypothesys
//- find preliminary besttrack as a reference
Float_t minchi2=10000;
Int_t maxn=0;
- AliITStrackV2 * besttrack=0;
+ AliITStrackMI * besttrack=0;
for (Int_t itrack=0;itrack<array->GetEntriesFast();itrack++){
- AliITStrackV2 * track = (AliITStrackV2*)array->At(itrack);
+ AliITStrackMI * track = (AliITStrackMI*)array->At(itrack);
if (!track) continue;
Float_t chi2 = NormalizedChi2(track,0);
//
}
}
else{
- delete array->RemoveAt(itrack);
+ if (track->fConstrain || track->fN>5){ //keep best short tracks - without vertex constrain
+ delete array->RemoveAt(itrack);
+ }
}
}
if (!besttrack) return;
Int_t * index = new Int_t[entries];
for (Int_t i=0;i<entries;i++) chi2[i] =10000;
for (Int_t itrack=0;itrack<entries;itrack++){
- AliITStrackV2 * track = (AliITStrackV2*)array->At(itrack);
+ AliITStrackMI * track = (AliITStrackMI*)array->At(itrack);
if (track){
track->fChi2MIP[0] = GetNormalizedChi2(track, mode);
if (track->fChi2MIP[0]<kMaxChi2PerCluster[0])
chi2[itrack] = track->fChi2MIP[0];
- else
- delete array->RemoveAt(itrack);
+ else{
+ if (track->fConstrain || track->fN>5){ //keep best short tracks - without vertex constrain
+ delete array->RemoveAt(itrack);
+ }
+ }
}
}
//
TMath::Sort(entries,chi2,index,kFALSE);
- besttrack = (AliITStrackV2*)array->At(index[0]);
+ besttrack = (AliITStrackMI*)array->At(index[0]);
if (besttrack&&besttrack->fChi2MIP[0]<kMaxChi2PerCluster[0]){
for (Int_t i=0;i<6;i++){
if (besttrack->fClIndex[i]>0){
// calculate one more time with updated normalized errors
for (Int_t i=0;i<entries;i++) chi2[i] =10000;
for (Int_t itrack=0;itrack<entries;itrack++){
- AliITStrackV2 * track = (AliITStrackV2*)array->At(itrack);
+ AliITStrackMI * track = (AliITStrackMI*)array->At(itrack);
if (track){
track->fChi2MIP[0] = GetNormalizedChi2(track,mode);
if (track->fChi2MIP[0]<kMaxChi2PerCluster[0])
chi2[itrack] = track->fChi2MIP[0]-0*(track->GetNumberOfClusters()+track->fNDeadZone);
- else
- delete array->RemoveAt(itrack);
+ else
+ {
+ if (track->fConstrain || track->fN>5){ //keep best short tracks - without vertex constrain
+ delete array->RemoveAt(itrack);
+ }
+ }
}
}
entries = array->GetEntriesFast();
//
+ //
if (entries>0){
TObjArray * newarray = new TObjArray();
TMath::Sort(entries,chi2,index,kFALSE);
- besttrack = (AliITStrackV2*)array->At(index[0]);
+ besttrack = (AliITStrackMI*)array->At(index[0]);
if (besttrack){
//
for (Int_t i=0;i<6;i++){
Float_t minn = besttrack->GetNumberOfClusters()-3;
Int_t accepted=0;
for (Int_t i=0;i<entries;i++){
- AliITStrackV2 * track = (AliITStrackV2*)array->At(index[i]);
+ AliITStrackMI * track = (AliITStrackMI*)array->At(index[i]);
if (!track) continue;
if (accepted>maxcut) break;
track->fChi2MIP[0] = GetNormalizedChi2(track,mode);
- if (track->GetNumberOfClusters()<6 && (track->fChi2MIP[0]+track->fNUsed>minchi2)){
- delete array->RemoveAt(index[i]);
- continue;
+ if (track->fConstrain || track->fN>5){ //keep best short tracks - without vertex constrain
+ if (track->GetNumberOfClusters()<6 && (track->fChi2MIP[0]+track->fNUsed>minchi2)){
+ delete array->RemoveAt(index[i]);
+ continue;
+ }
}
- if (track->fChi2MIP[0]+track->fNUsed<minchi2 && track->GetNumberOfClusters()>=minn){
- accepted++;
+ Bool_t shortbest = !track->fConstrain && track->fN<6;
+ if ((track->fChi2MIP[0]+track->fNUsed<minchi2 && track->GetNumberOfClusters()>=minn) ||shortbest){
+ if (!shortbest) accepted++;
//
newarray->AddLast(array->RemoveAt(index[i]));
for (Int_t i=0;i<6;i++){
-AliITStrackV2 * AliITStrackerMI::GetBestHypothesys(Int_t esdindex, AliITStrackV2 * original, Int_t checkmax)
+AliITStrackMI * AliITStrackerMI::GetBestHypothesys(Int_t esdindex, AliITStrackMI * original, Int_t checkmax)
{
//-------------------------------------------------------------
// try to find best hypothesy
Int_t entries = array->GetEntriesFast();
if (!entries) return 0;
Float_t minchi2 = 100000;
- AliITStrackV2 * besttrack=0;
+ AliITStrackMI * besttrack=0;
//
- AliITStrackV2 * backtrack = new AliITStrackV2(*original);
- AliITStrackV2 * forwardtrack = new AliITStrackV2(*original);
+ AliITStrackMI * backtrack = new AliITStrackMI(*original);
+ AliITStrackMI * forwardtrack = new AliITStrackMI(*original);
+ Double_t xyzv[]={GetX(),GetY(),GetZ()};
+ Double_t ersv[]={GetSigmaX()/3.,GetSigmaY()/3.,GetSigmaZ()/3.};
//
for (Int_t i=0;i<entries;i++){
- AliITStrackV2 * track = (AliITStrackV2*)array->At(i);
+ AliITStrackMI * track = (AliITStrackMI*)array->At(i);
if (!track) continue;
+ Float_t sigmarfi,sigmaz;
+ GetDCASigma(track,sigmarfi,sigmaz);
+ track->fDnorm[0] = sigmarfi;
+ track->fDnorm[1] = sigmaz;
+ //
track->fChi2MIP[1] = 1000000;
track->fChi2MIP[2] = 1000000;
track->fChi2MIP[3] = 1000000;
//
// backtrack
- backtrack = new(backtrack) AliITStrackV2(*track);
- backtrack->ResetCovariance();
- backtrack->ResetCovariance();
+ backtrack = new(backtrack) AliITStrackMI(*track);
+ if (track->fConstrain){
+ if (!backtrack->PropagateTo(3.,0.0028,65.19)) continue;
+ if (!backtrack->Improve(0,xyzv,ersv)) continue;
+ if (!backtrack->PropagateTo(2.,0.0028,0)) continue;
+ if (!backtrack->Improve(0,xyzv,ersv)) continue;
+ if (!backtrack->PropagateTo(1.,0.0028,0)) continue;
+ if (!backtrack->Improve(0,xyzv,ersv)) continue;
+ if (!backtrack->PropagateToVertex()) continue;
+ backtrack->ResetCovariance();
+ if (!backtrack->Improve(0,xyzv,ersv)) continue;
+ }else{
+ backtrack->ResetCovariance();
+ }
backtrack->ResetClusters();
+
Double_t x = original->GetX();
if (!RefitAt(x,backtrack,track)) continue;
+ //
track->fChi2MIP[1] = NormalizedChi2(backtrack,0);
//for (Int_t i=2;i<6;i++){track->fDy[i]+=backtrack->fDy[i]; track->fDz[i]+=backtrack->fDz[i];}
if (track->fChi2MIP[1]>kMaxChi2PerCluster[1]*6.) continue;
if (!(track->fConstrain)&&track->fChi2MIP[1]>kMaxChi2PerCluster[1]) continue;
Bool_t isOK=kTRUE;
- /*
- for (Int_t i=0;i<6;i++){
- if (track->fClIndex[i]>0){
- Double_t dy1 = (track->fDy[i]/track->fSigmaY[i]);
- Double_t dz1 = (track->fDz[i]/track->fSigmaZ[i]);
- Double_t dy2 = (backtrack->fDy[i]/backtrack->fSigmaY[i]);
- Double_t dz2 = (backtrack->fDz[i]/backtrack->fSigmaZ[i]);
- if (TMath::Min(dy1*dy1+dz1*dz1,dy2*dy2+dz2*dz2)> kMaxChi2sR[i]) isOK =kFALSE;
- track->fDy[i+6] = backtrack->fDy[i];
- track->fDz[i+6] = backtrack->fDz[i];
- track->fSigmaY[i+6] = backtrack->fSigmaY[i];
- track->fSigmaZ[i+6] = backtrack->fSigmaZ[i];
- }
- else{
- if (i==5){
- if (track->fClIndex[i-1]>0){
- Double_t dy2 = (backtrack->fDy[i-1]/backtrack->fSigmaY[i-1]);
- Double_t dz2 = (backtrack->fDz[i-1]/backtrack->fSigmaZ[i-1]);
- if (dy2*dy2+dz2*dz2> kMaxChi2sR[i]) isOK =kFALSE;
- }
- else isOK = kFALSE;
- }
- }
- }
- */
if(!isOK) continue;
//
//forward track - without constraint
- forwardtrack = new(forwardtrack) AliITStrackV2(*original);
+ forwardtrack = new(forwardtrack) AliITStrackMI(*original);
forwardtrack->ResetClusters();
x = track->GetX();
- if (!RefitAt(x,forwardtrack,track)) continue;
+ RefitAt(x,forwardtrack,track);
track->fChi2MIP[2] = NormalizedChi2(forwardtrack,0);
if (track->fChi2MIP[2]>kMaxChi2PerCluster[2]*6.0) continue;
if (!(track->fConstrain)&&track->fChi2MIP[2]>kMaxChi2PerCluster[2]) continue;
forwardtrack->fChi2MIP[ichi] = track->fChi2MIP[ichi];
}
if (chi2 < minchi2){
- //besttrack = new AliITStrackV2(*forwardtrack);
+ //besttrack = new AliITStrackMI(*forwardtrack);
besttrack = track;
besttrack->SetLabel(track->GetLabel());
besttrack->fFakeRatio = track->fFakeRatio;
delete forwardtrack;
Int_t accepted=0;
for (Int_t i=0;i<entries;i++){
- AliITStrackV2 * track = (AliITStrackV2*)array->At(i);
+ AliITStrackMI * track = (AliITStrackMI*)array->At(i);
if (!track) continue;
+
if (accepted>checkmax || track->fChi2MIP[3]>kMaxChi2PerCluster[3]*6. ||
(track->GetNumberOfClusters()<besttrack->GetNumberOfClusters()-1.)||
track->fChi2MIP[0]>besttrack->fChi2MIP[0]+2.*besttrack->fNUsed+3.){
- delete array->RemoveAt(i);
- continue;
+ if (track->fConstrain || track->fN>5){ //keep best short tracks - without vertex constrain
+ delete array->RemoveAt(i);
+ continue;
+ }
}
else{
accepted++;
array->Compress();
SortTrackHypothesys(esdindex,checkmax,1);
array = (TObjArray*) fTrackHypothesys.At(esdindex);
- besttrack = (AliITStrackV2*)array->At(0);
+ besttrack = (AliITStrackMI*)array->At(0);
if (!besttrack) return 0;
besttrack->fChi2MIP[8]=0;
fBestTrackIndex[esdindex]=0;
entries = array->GetEntriesFast();
- AliITStrackV2 *longtrack =0;
+ AliITStrackMI *longtrack =0;
minchi2 =1000;
Float_t minn=besttrack->GetNumberOfClusters()+besttrack->fNDeadZone;
for (Int_t itrack=entries-1;itrack>0;itrack--){
- AliITStrackV2 * track = (AliITStrackV2*)array->At(itrack);
+ AliITStrackMI * track = (AliITStrackMI*)array->At(itrack);
if (!track->fConstrain) continue;
if (track->GetNumberOfClusters()+track->fNDeadZone<minn) continue;
if (track->fChi2MIP[0]-besttrack->fChi2MIP[0]>0.0) continue;
Float_t shared = GetNumberOfSharedClusters(besttrack,esdindex,list,clist);
if (besttrack->fConstrain&&besttrack->fChi2MIP[0]<kMaxChi2PerCluster[0]&&besttrack->fChi2MIP[1]<kMaxChi2PerCluster[1]
&&besttrack->fChi2MIP[2]<kMaxChi2PerCluster[2]&&besttrack->fChi2MIP[3]<kMaxChi2PerCluster[3]){
- //if (besttrack->fChi2MIP[0]<kMaxChi2PerCluster[0]&&besttrack->fChi2MIP[1]<kMaxChi2PerCluster[1]
- // &&besttrack->fChi2MIP[2]<kMaxChi2PerCluster[2]&&besttrack->fChi2MIP[3]<kMaxChi2PerCluster[3]){
RegisterClusterTracks(besttrack,esdindex);
}
//
if (shared>2.5) return 0;
if (shared>1.0) return besttrack;
//
- // don't sign clusters if not gold track
- Double_t deltad = besttrack->GetD(GetX(),GetY());
- Double_t deltaz = besttrack->GetZat(GetX()) - GetZ();
- Double_t deltaprim = TMath::Sqrt(deltad*deltad+deltaz*deltaz);
- if (deltaprim>0.1 && (fConstraint[fPass])) return besttrack;
- if (TMath::Abs(deltad)>0.1) return besttrack;
- if (TMath::Abs(deltaz)>0.1) return besttrack;
- if (besttrack->fChi2MIP[0]>4.) return besttrack;
- if (besttrack->fChi2MIP[1]>4.) return besttrack;
- if (besttrack->fChi2MIP[2]>4.) return besttrack;
- if (besttrack->fChi2MIP[3]>4.) return besttrack;
+ // Don't sign clusters if not gold track
+ //
+ if (!besttrack->IsGoldPrimary()) return besttrack;
+ if (besttrack->fESDtrack->GetKinkIndex(0)!=0) return besttrack; //track belong to kink
+ //
if (fConstraint[fPass]){
//
// sign clusters
Bool_t cansign = kTRUE;
for (Int_t itrack=0;itrack<entries; itrack++){
- AliITStrackV2 * track = (AliITStrackV2*)array->At(i);
+ AliITStrackMI * track = (AliITStrackMI*)array->At(i);
if (!track) continue;
if (track->fChi2MIP[0]>besttrack->fChi2MIP[0]+2.*shared+1.) break;
if ( (track->fClIndex[ilayer]>0) && (track->fClIndex[ilayer]!=besttrack->fClIndex[ilayer])){
{
//
// get "best" hypothesys
- //for (Int_t ilayer=0;ilayer<6;ilayer++) fgLayers[ilayer].ResetWeights();
-
+ //
Int_t nentries = itsTracks.GetEntriesFast();
for (Int_t i=0;i<nentries;i++){
- AliITStrackV2* track = (AliITStrackV2*)itsTracks.At(i);
+ AliITStrackMI* track = (AliITStrackMI*)itsTracks.At(i);
if (!track) continue;
TObjArray * array = (TObjArray*) fTrackHypothesys.At(i);
if (!array) continue;
if (array->GetEntriesFast()<=0) continue;
//
- AliITStrackV2* longtrack=0;
+ AliITStrackMI* longtrack=0;
Float_t minn=0;
Float_t maxchi2=1000;
for (Int_t j=0;j<array->GetEntriesFast();j++){
- AliITStrackV2* track = (AliITStrackV2*)array->At(j);
+ AliITStrackMI* track = (AliITStrackMI*)array->At(j);
if (!track) continue;
+ if (track->fGoldV0) {
+ longtrack = track; //gold V0 track taken
+ break;
+ }
if (track->GetNumberOfClusters()+track->fNDeadZone<minn) continue;
- if (track->GetNumberOfClusters()+track->fNDeadZone>minn) maxchi2 = track->fChi2MIP[0];
- if (track->fChi2MIP[0]>maxchi2) continue;
+ Float_t chi2 = track->fChi2MIP[0];
+ if (fAfterV0){
+ if (!track->fGoldV0&&track->fConstrain==kFALSE) chi2+=5;
+ }
+ if (track->GetNumberOfClusters()+track->fNDeadZone>minn) maxchi2 = track->fChi2MIP[0];
+ //
+ if (chi2 > maxchi2) continue;
minn= track->GetNumberOfClusters()+track->fNDeadZone;
- maxchi2 = track->fChi2MIP[0];
+ maxchi2 = chi2;
longtrack=track;
- break;
}
- AliITStrackV2 * besttrack = (AliITStrackV2*)array->At(0);
+ //
+ //
+ //
+ AliITStrackMI * besttrack = (AliITStrackMI*)array->At(0);
if (!longtrack) {longtrack = besttrack;}
else besttrack= longtrack;
+ //
if (besttrack){
Int_t list[6];
AliITSclusterV2 * clist[6];
besttrack = GetBest2Tracks(i,sharedtrack,10,5.5);
}
}
+ if (besttrack&&fAfterV0){
+ UpdateESDtrack(besttrack,AliESDtrack::kITSin);
+ }
if (besttrack&&fConstraint[fPass])
UpdateESDtrack(besttrack,AliESDtrack::kITSin);
//if (besttrack&&besttrack->fConstrain)
}
-void AliITStrackerMI::CookLabel(AliITStrackV2 *track,Float_t wrong) const {
+void AliITStrackerMI::CookLabel(AliITStrackMI *track,Float_t wrong) const {
//--------------------------------------------------------------------
//This function "cooks" a track label. If label<0, this track is fake.
//--------------------------------------------------------------------
-void AliITStrackerMI::CookdEdx(AliITStrackV2* track)
+void AliITStrackerMI::CookdEdx(AliITStrackMI* track)
{
//
//
void AliITStrackerMI::MakeCoeficients(Int_t ntracks){
//
//
+ if (fCoeficients) delete []fCoeficients;
fCoeficients = new Float_t[ntracks*48];
for (Int_t i=0;i<ntracks*48;i++) fCoeficients[i]=-1.;
}
-Double_t AliITStrackerMI::GetPredictedChi2MI(AliITStrackV2* track, const AliITSclusterV2 *cluster,Int_t layer)
+Double_t AliITStrackerMI::GetPredictedChi2MI(AliITStrackMI* track, const AliITSclusterV2 *cluster,Int_t layer)
{
//
//
}
-Int_t AliITStrackerMI::UpdateMI(AliITStrackV2* track, const AliITSclusterV2* cl,Double_t chi2,Int_t index) const
+Int_t AliITStrackerMI::UpdateMI(AliITStrackMI* track, const AliITSclusterV2* cl,Double_t chi2,Int_t index) const
{
//
//
//STRIPS
if (layer>3){
+ //factor 1.8 appears in new simulation
+ //
+ Float_t scale=1.8;
if (cl->GetNy()==100||cl->GetNz()==100){
- erry = 0.004;
- errz = 0.2;
+ erry = 0.004*scale;
+ errz = 0.2*scale;
return 100;
}
if (cl->GetNy()+cl->GetNz()>12){
- erry = 0.06;
- errz = 0.57;
+ erry = 0.06*scale;
+ errz = 0.57*scale;
return 100;
}
Float_t normq = cl->GetQ()/(TMath::Sqrt(1+theta*theta+phi*phi));
//
if (cl->GetType()==1 || cl->GetType()==10 ){
if (chargematch<1.0 || (cl->GetNy()+cl->GetNz()<nz+ny+0.5)){
- errz = 0.043;
- erry = 0.00094;
+ errz = 0.043*scale;
+ erry = 0.00094*scale;
return 101;
}
if (cl->GetNy()+cl->GetNz()<nz+ny+1.2){
- errz = 0.06;
- erry =0.0013;
+ errz = 0.06*scale;
+ erry =0.0013*scale;
return 102;
}
- erry = 0.0027;
- errz = TMath::Min(0.028*(chargematch+cl->GetNy()+cl->GetNz()-nz+ny),0.15);
+ erry = 0.0027*scale;
+ errz = TMath::Min(0.028*(chargematch+cl->GetNy()+cl->GetNz()-nz+ny),0.15)*scale;
return 103;
}
if (cl->GetType()==2 || cl->GetType()==11 ){
- erry = TMath::Min(0.0010*(1+chargematch+cl->GetNy()+cl->GetNz()-nz+ny),0.05);
- errz = TMath::Min(0.025*(1+chargematch+cl->GetNy()+cl->GetNz()-nz+ny),0.5);
+ erry = TMath::Min(0.0010*(1+chargematch+cl->GetNy()+cl->GetNz()-nz+ny),0.05)*scale;
+ errz = TMath::Min(0.025*(1+chargematch+cl->GetNy()+cl->GetNz()-nz+ny),0.5)*scale;
return 104;
}
if (cl->GetType()>100 ){
if ((chargematch+cl->GetNy()+cl->GetNz()-nz-ny<1.5)){
- errz = 0.05;
- erry = 0.00096;
+ errz = 0.05*scale;
+ erry = 0.00096*scale;
return 105;
}
if (cl->GetNy()+cl->GetNz()-nz-ny<1){
- errz = 0.10;
- erry = 0.0025;
+ errz = 0.10*scale;
+ erry = 0.0025*scale;
return 106;
}
- errz = TMath::Min(0.05*(chargematch+cl->GetNy()+cl->GetNz()-nz-ny),0.4);
- erry = TMath::Min(0.003*(chargematch+cl->GetNy()+cl->GetNz()-nz-ny),0.05);
+ errz = TMath::Min(0.05*(chargematch+cl->GetNy()+cl->GetNz()-nz-ny),0.4)*scale;
+ erry = TMath::Min(0.003*(chargematch+cl->GetNy()+cl->GetNz()-nz-ny),0.05)*scale;
return 107;
}
Float_t diff = cl->GetNy()+cl->GetNz()-ny-nz;
+void AliITStrackerMI::GetDCASigma(AliITStrackMI* track, Float_t & sigmarfi, Float_t &sigmaz)
+{
+ //
+ //DCA sigmas parameterization
+ //to be paramterized using external parameters in future
+ //
+ //
+ sigmarfi = 0.004+1.4 *TMath::Abs(track->fP4)+332.*track->fP4*track->fP4;
+ sigmaz = 0.011+4.37*TMath::Abs(track->fP4);
+}
+
+
void AliITStrackerMI::SignDeltas( TObjArray *ClusterArray, Float_t vz)
{
//
}
-void AliITStrackerMI::UpdateESDtrack(AliITStrackV2* track, ULong_t flags) const
+void AliITStrackerMI::UpdateESDtrack(AliITStrackMI* track, ULong_t flags) const
{
//
//
track->UpdateESDtrack(flags);
- AliITStrackV2 * oldtrack = (AliITStrackV2*)(track->fESDtrack->GetITStrack());
+ AliITStrackMI * oldtrack = (AliITStrackMI*)(track->fESDtrack->GetITStrack());
if (oldtrack) delete oldtrack;
- track->fESDtrack->SetITStrack(new AliITStrackV2(*track));
+ track->fESDtrack->SetITStrack(new AliITStrackMI(*track));
+ if (TMath::Abs(track->fDnorm[1])<0.000000001){
+ printf("Problem\n");
+ }
}
-
-void AliITStrackerMI::FindV0(AliESD */*event*/)
-{
- //
- // fast V0 finder
+Int_t AliITStrackerMI::GetNearestLayer(const Double_t *xr) const{
//
- AliHelix helixes[30000];
- TObjArray trackarray(30000);
- Float_t dist[30000];
- AliITSRecV0Info vertex;
+ //Get nearest upper layer close to the point xr.
+ // rough approximation
//
- Int_t entries = fTrackHypothesys.GetEntriesFast();
- for (Int_t i=0;i<entries;i++){
- TObjArray * array = (TObjArray*)fTrackHypothesys.At(i);
- if (!array) continue;
- AliITStrackV2 * track = (AliITStrackV2*)array->At(fBestTrackIndex[i]);
- if (track){
- dist[i] = TMath::Sqrt(track->fD[0]*track->fD[0]+track->fD[1]*track->fD[1]);
- trackarray.AddAt(track,i);
- new (&helixes[i]) AliHelix(*track);
- }
- }
- for (Int_t itrack0=0;itrack0<entries;itrack0++){
- AliITStrackV2 * track0 = (AliITStrackV2*)trackarray.At(itrack0);
- if (!track0) continue;
- if (dist[itrack0]<0.2) continue;
- for (Int_t itrack1=itrack0+1;itrack1<entries;itrack1++){
- AliITStrackV2 * track1 = (AliITStrackV2*)trackarray.At(itrack1);
- if (!track1) continue;
- if (dist[itrack1]<0.2) continue;
- if (track1->fP4*track0->fP4>0) continue; //the same sign
- AliHelix *h1 = &helixes[itrack0];
- AliHelix *h2 = &helixes[itrack1];
-
- Double_t distance = TestV0(h1,h2,&vertex);
- if (distance>0.4) continue;
- if (vertex.fRr<0.3) continue;
- if (vertex.fRr>27) continue;
- Float_t v[3]={GetX(),GetY(),GetZ()};
- vertex.Update(v,track0->fMass, track1->fMass);
-
- if ((TMath::Abs(vertex.fInvMass-0.4976)<0.05 || TMath::Abs(vertex.fInvMass-0.0)<0.1|| TMath::Abs(vertex.fInvMass-1.1)<0.1)) {
- if (vertex.fPointAngle<0.85) continue;
- }
- else{
- if (vertex.fPointAngle<0.98) continue;
- }
- if (TMath::Abs(TMath::Abs(track0->fLab)-TMath::Abs(track1->fLab))<2){
- Float_t mindist = FindBestPair(itrack0,itrack1,&vertex);
- if (mindist>1) FindBestPair(itrack0,itrack1,&vertex);
- }
- }
- }
-}
-
-Double_t AliITStrackerMI::FindBestPair(Int_t esdtrack0, Int_t esdtrack1,AliITSRecV0Info *vertex)
-{
- //
- // try to find best pair from the tree of track hyp.
- //
- TObjArray * array0 = (TObjArray*)fTrackHypothesys.At(esdtrack0);
- Int_t entries0 = array0->GetEntriesFast();
- TObjArray * array1 = (TObjArray*)fTrackHypothesys.At(esdtrack1);
- Int_t entries1 = array1->GetEntriesFast();
- //
- //
- AliITSRecV0Info v0;
- AliITSRecV0Info vbest;
- Double_t criticalradius = vertex->fRr;
- Double_t mindist =1000;
- Int_t bestpair[2];
- //
- for (Int_t itrack0=0;itrack0<entries0;itrack0++){
- AliITStrackV2 * track0 = (AliITStrackV2*)array0->At(itrack0);
- if (!track0) continue;
- if (track0->fX<criticalradius-1) continue;
- if (track0->fX>criticalradius+5) continue;
- for (Int_t itrack1=0;itrack1<entries1;itrack1++){
- AliITStrackV2 * track1 = (AliITStrackV2*)array1->At(itrack1);
- if (!track1) continue;
- if (track1->fX<criticalradius-1) continue;
- if (track1->fX>criticalradius+5) continue;
-
- AliHelix h0(*track0);
- AliHelix h1(*track1);
- Double_t distance = TestV0(&h0,&h1,&v0);
- if (v0.fRr>30) continue;
- if (v0.fRr<0.2) continue;
- Float_t v[3]={GetX(),GetY(),GetZ()};
- v0.Update(v,track0->fMass, track1->fMass);
- if (distance<mindist){
- mindist = distance;
- bestpair[0]=itrack0;
- bestpair[1]=itrack1;
- vbest = v0;
- }
+ const Float_t radiuses[6]={4,6.5,15.03,24.,38.5,43.7};
+ Float_t radius = TMath::Sqrt(xr[0]*xr[0]+xr[1]*xr[1]);
+ Int_t res =6;
+ for (Int_t i=0;i<6;i++){
+ if (radius<radiuses[i]){
+ res =i;
+ break;
}
}
- if (mindist<0.2){
- vbest.Dump();
- vertex->Dump();
- }
- return mindist;
+ return res;
}
-void AliITSRecV0Info::Update(Float_t vertex[3], Float_t mass1, Float_t mass2)
-{
+void AliITStrackerMI::UpdateTPCV0(AliESD *event){
//
- //Update V0 information
+ //try to update, or reject TPC V0s
//
- Float_t v[3] = {fXr[0]-vertex[0],fXr[1]-vertex[1],fXr[2]-vertex[2]};
- Float_t p[3] = {fPdr[0]+fPm[0], fPdr[1]+fPm[1],fPdr[2]+fPm[2]};
-
- Float_t vnorm2 = v[0]*v[0]+v[1]*v[1];
- Float_t vnorm3 = TMath::Sqrt(v[2]*v[2]+vnorm2);
- vnorm2 = TMath::Sqrt(vnorm2);
- Float_t pnorm2 = p[0]*p[0]+p[1]*p[1];
- Float_t pnorm3 = TMath::Sqrt(p[2]*p[2]+pnorm2);
- pnorm2 = TMath::Sqrt(pnorm2);
-
- fPointAngleFi = (v[0]*p[0]+v[1]*p[1])/(vnorm2*pnorm2);
- fPointAngleTh = (v[2]*p[2]+vnorm2*pnorm2)/(vnorm3*pnorm3);
- fPointAngle = (v[0]*p[0]+v[1]*p[1]+v[2]*p[2])/(vnorm3*pnorm3);
-
- Float_t e1 = TMath::Sqrt(mass1*mass1+
- fPdr[0]*fPdr[0]+
- fPdr[1]*fPdr[1]+
- fPdr[2]*fPdr[2]);
- Float_t e2 = TMath::Sqrt(mass2*mass2+
- fPm[0]*fPm[0]+
- fPm[1]*fPm[1]+
- fPm[2]*fPm[2]);
-
- fInvMass =
- (fPm[0]+fPdr[0])*(fPm[0]+fPdr[0])+
- (fPm[1]+fPdr[1])*(fPm[1]+fPdr[1])+
- (fPm[2]+fPdr[2])*(fPm[2]+fPdr[2]);
-
- fInvMass = TMath::Sqrt((e1+e2)*(e1+e2)-fInvMass);
+ Int_t nv0s = event->GetNumberOfV0MIs();
+ Int_t nitstracks = fTrackHypothesys.GetEntriesFast();
+ for (Int_t i=0;i<nv0s;i++){
+ AliESDV0MI * vertex = event->GetV0MI(i);
+ Int_t ip = vertex->GetIndex(0);
+ Int_t im = vertex->GetIndex(1);
+ //
+ TObjArray * arrayp = (ip<nitstracks) ? (TObjArray*)fTrackHypothesys.At(ip):0;
+ TObjArray * arraym = (im<nitstracks) ? (TObjArray*)fTrackHypothesys.At(im):0;
+ AliITStrackMI * trackp = (arrayp!=0) ? (AliITStrackMI*)arrayp->At(0):0;
+ AliITStrackMI * trackm = (arraym!=0) ? (AliITStrackMI*)arraym->At(0):0;
+ //
+ //
+ if (trackp){
+ if (trackp->fN+trackp->fNDeadZone>5.5){
+ if (trackp->fConstrain&&trackp->fChi2MIP[0]<3) vertex->SetStatus(-100);
+ if (!trackp->fConstrain&&trackp->fChi2MIP[0]<2) vertex->SetStatus(-100);
+ }
+ }
+ if (trackm){
+ if (trackm->fN+trackm->fNDeadZone>5.5){
+ if (trackm->fConstrain&&trackm->fChi2MIP[0]<3) vertex->SetStatus(-100);
+ if (!trackm->fConstrain&&trackm->fChi2MIP[0]<2) vertex->SetStatus(-100);
+ }
+ }
+ if (vertex->GetStatus()==-100) continue;
+ //
+ Int_t clayer = GetNearestLayer(vertex->GetXrp());
+ vertex->SetNBefore(clayer); //
+ vertex->SetChi2Before(9*clayer); //
+ vertex->SetNAfter(6-clayer); //
+ vertex->SetChi2After(0); //
+ //
+ if (clayer >1 ){ // calculate chi2 before vertex
+ Float_t chi2p = 0, chi2m=0;
+ //
+ if (trackp){
+ for (Int_t ilayer=0;ilayer<clayer;ilayer++){
+ if (trackp->fClIndex[ilayer]>0){
+ chi2p+=trackp->fDy[ilayer]*trackp->fDy[ilayer]/(trackp->fSigmaY[ilayer]*trackp->fSigmaY[ilayer])+
+ trackp->fDz[ilayer]*trackp->fDz[ilayer]/(trackp->fSigmaZ[ilayer]*trackp->fSigmaZ[ilayer]);
+ }
+ else{
+ chi2p+=9;
+ }
+ }
+ }else{
+ chi2p = 9*clayer;
+ }
+ //
+ if (trackm){
+ for (Int_t ilayer=0;ilayer<clayer;ilayer++){
+ if (trackm->fClIndex[ilayer]>0){
+ chi2m+=trackm->fDy[ilayer]*trackm->fDy[ilayer]/(trackm->fSigmaY[ilayer]*trackm->fSigmaY[ilayer])+
+ trackm->fDz[ilayer]*trackm->fDz[ilayer]/(trackm->fSigmaZ[ilayer]*trackm->fSigmaZ[ilayer]);
+ }
+ else{
+ chi2m+=9;
+ }
+ }
+ }else{
+ chi2m = 9*clayer;
+ }
+ vertex->SetChi2Before(TMath::Min(chi2p,chi2m));
+ if (TMath::Min(chi2p,chi2m)/Float_t(clayer)<4) vertex->SetStatus(-10); // track exist before vertex
+ }
+
+ if (clayer < 5 ){ // calculate chi2 after vertex
+ Float_t chi2p = 0, chi2m=0;
+ //
+ if (trackp&&TMath::Abs(trackp->fP3)<1.){
+ for (Int_t ilayer=clayer;ilayer<6;ilayer++){
+ if (trackp->fClIndex[ilayer]>0){
+ chi2p+=trackp->fDy[ilayer]*trackp->fDy[ilayer]/(trackp->fSigmaY[ilayer]*trackp->fSigmaY[ilayer])+
+ trackp->fDz[ilayer]*trackp->fDz[ilayer]/(trackp->fSigmaZ[ilayer]*trackp->fSigmaZ[ilayer]);
+ }
+ else{
+ chi2p+=9;
+ }
+ }
+ }else{
+ chi2p = 0;
+ }
+ //
+ if (trackm&&TMath::Abs(trackm->fP3)<1.){
+ for (Int_t ilayer=clayer;ilayer<6;ilayer++){
+ if (trackm->fClIndex[ilayer]>0){
+ chi2m+=trackm->fDy[ilayer]*trackm->fDy[ilayer]/(trackm->fSigmaY[ilayer]*trackm->fSigmaY[ilayer])+
+ trackm->fDz[ilayer]*trackm->fDz[ilayer]/(trackm->fSigmaZ[ilayer]*trackm->fSigmaZ[ilayer]);
+ }
+ else{
+ chi2m+=9;
+ }
+ }
+ }else{
+ chi2m = 0;
+ }
+ vertex->SetChi2After(TMath::Max(chi2p,chi2m));
+ if (TMath::Max(chi2m,chi2p)/Float_t(6-clayer)>9) vertex->SetStatus(-20); // track not found in ITS
+ }
+ }
+ //
}
-
-
-Double_t AliITStrackerMI::TestV0(AliHelix *helix1, AliHelix *helix2, AliITSRecV0Info *vertex)
+void AliITStrackerMI::FindV02(AliESD *event)
{
//
- // test the helixes for the distnce calculate vertex characteristic
- //
- Float_t distance1,distance2;
- AliHelix & dhelix1 = *helix1;
- Double_t pp[3],xx[3];
- dhelix1.GetMomentum(0,pp,0);
- dhelix1.Evaluate(0,xx);
- AliHelix &mhelix = *helix2;
+ // V0 finder
+ //
+ // Cuts on DCA - R dependent
+ // max distance DCA between 2 tracks cut
+ // maxDist = TMath::Min(kMaxDist,kMaxDist0+pvertex->GetRr()*kMaxDist);
+ //
+ const Float_t kMaxDist0 = 0.1;
+ const Float_t kMaxDist1 = 0.1;
+ const Float_t kMaxDist = 1;
+ const Float_t kMinPointAngle = 0.85;
+ const Float_t kMinPointAngle2 = 0.99;
+ const Float_t kMinR = 0.5;
+ const Float_t kMaxR = 220;
+ //const Float_t kCausality0Cut = 0.19;
+ //const Float_t kLikelihood01Cut = 0.25;
+ //const Float_t kPointAngleCut = 0.9996;
+ const Float_t kCausality0Cut = 0.19;
+ const Float_t kLikelihood01Cut = 0.45;
+ const Float_t kLikelihood1Cut = 0.5;
+ const Float_t kCombinedCut = 0.55;
+
+ //
+ //
+ TTreeSRedirector &cstream = *fDebugStreamer;
+ Int_t ntracks = event->GetNumberOfTracks();
+ Int_t nitstracks = fTrackHypothesys.GetEntriesFast();
+ fOriginal.Expand(ntracks);
+ fTrackHypothesys.Expand(ntracks);
+ fBestHypothesys.Expand(ntracks);
+ //
+ AliHelix * helixes = new AliHelix[ntracks+2];
+ TObjArray trackarray(ntracks+2); //array with tracks - with vertex constrain
+ TObjArray trackarrayc(ntracks+2); //array of "best tracks" - without vertex constrain
+ TObjArray trackarrayl(ntracks+2); //array of "longest tracks" - without vertex constrain
+ Bool_t * forbidden = new Bool_t [ntracks+2];
+ Int_t *itsmap = new Int_t [ntracks+2];
+ Float_t *dist = new Float_t[ntracks+2];
+ Float_t *normdist0 = new Float_t[ntracks+2];
+ Float_t *normdist1 = new Float_t[ntracks+2];
+ Float_t *normdist = new Float_t[ntracks+2];
+ Float_t *norm = new Float_t[ntracks+2];
+ Float_t *maxr = new Float_t[ntracks+2];
+ Float_t *minr = new Float_t[ntracks+2];
+ Float_t *minPointAngle= new Float_t[ntracks+2];
+ //
+ AliESDV0MI *pvertex = new AliESDV0MI;
+ AliITStrackMI * dummy= new AliITStrackMI;
+ dummy->SetLabel(0);
+ AliITStrackMI trackat0; //temporary track for DCA calculation
+ //
+ Float_t primvertex[3]={GetX(),GetY(),GetZ()};
+ //
+ // make its - esd map
+ //
+ for (Int_t itrack=0;itrack<ntracks+2;itrack++) {
+ itsmap[itrack] = -1;
+ forbidden[itrack] = kFALSE;
+ maxr[itrack] = kMaxR;
+ minr[itrack] = kMinR;
+ minPointAngle[itrack] = kMinPointAngle;
+ }
+ for (Int_t itrack=0;itrack<nitstracks;itrack++){
+ AliITStrackMI * original = (AliITStrackMI*)(fOriginal.At(itrack));
+ Int_t esdindex = original->fESDtrack->GetID();
+ itsmap[esdindex] = itrack;
+ }
+ //
+ // create its tracks from esd tracks if not done before
+ //
+ for (Int_t itrack=0;itrack<ntracks;itrack++){
+ if (itsmap[itrack]>=0) continue;
+ AliITStrackMI * tpctrack = new AliITStrackMI(*(event->GetTrack(itrack)));
+ tpctrack->fD[0] = tpctrack->GetD(GetX(),GetY());
+ tpctrack->fD[1] = tpctrack->GetZat(GetX())-GetZ();
+ if (tpctrack->fD[0]<20 && tpctrack->fD[1]<20){
+ // tracks which can reach inner part of ITS
+ // propagate track to outer its volume - with correction for material
+ CorrectForDeadZoneMaterial(tpctrack);
+ }
+ itsmap[itrack] = nitstracks;
+ fOriginal.AddAt(tpctrack,nitstracks);
+ nitstracks++;
+ }
+ //
+ // fill temporary arrays
+ //
+ for (Int_t itrack=0;itrack<ntracks;itrack++){
+ AliESDtrack * esdtrack = event->GetTrack(itrack);
+ Int_t itsindex = itsmap[itrack];
+ AliITStrackMI *original = (AliITStrackMI*)fOriginal.At(itsmap[itrack]);
+ if (!original) continue;
+ AliITStrackMI *bestConst = 0;
+ AliITStrackMI *bestLong = 0;
+ AliITStrackMI *best = 0;
+ //
+ //
+ TObjArray * array = (TObjArray*) fTrackHypothesys.At(itsindex);
+ Int_t hentries = (array==0) ? 0 : array->GetEntriesFast();
+ // Get best track with vertex constrain
+ for (Int_t ih=0;ih<hentries;ih++){
+ AliITStrackMI * trackh = (AliITStrackMI*)array->At(ih);
+ if (!trackh->fConstrain) continue;
+ if (!bestConst) bestConst = trackh;
+ if (trackh->fN>5.0){
+ bestConst = trackh; // full track - with minimal chi2
+ break;
+ }
+ if (trackh->fN+trackh->fNDeadZone<=bestConst->fN+bestConst->fNDeadZone) continue;
+ bestConst = trackh;
+ break;
+ }
+ // Get best long track without vertex constrain and best track without vertex constrain
+ for (Int_t ih=0;ih<hentries;ih++){
+ AliITStrackMI * trackh = (AliITStrackMI*)array->At(ih);
+ if (trackh->fConstrain) continue;
+ if (!best) best = trackh;
+ if (!bestLong) bestLong = trackh;
+ if (trackh->fN>5.0){
+ bestLong = trackh; // full track - with minimal chi2
+ break;
+ }
+ if (trackh->fN+trackh->fNDeadZone<=bestLong->fN+bestLong->fNDeadZone) continue;
+ bestLong = trackh;
+ }
+ if (!best) {
+ best = original;
+ bestLong = original;
+ }
+ trackat0 = *bestLong;
+ Double_t xx,yy,zz,alpha;
+ bestLong->GetGlobalXYZat(bestLong->GetX(),xx,yy,zz);
+ alpha = TMath::ATan2(yy,xx);
+ trackat0.Propagate(alpha,0);
+ // calculate normalized distances to the vertex
+ //
+ if ( bestLong->fN>3 ){
+ dist[itsindex] = trackat0.fP0;
+ norm[itsindex] = TMath::Sqrt(trackat0.fC00);
+ normdist0[itsindex] = TMath::Abs(trackat0.fP0/norm[itsindex]);
+ normdist1[itsindex] = TMath::Abs((trackat0.fP1-primvertex[2])/TMath::Sqrt(trackat0.fC11));
+ normdist[itsindex] = TMath::Sqrt(normdist0[itsindex]*normdist0[itsindex]+normdist1[itsindex]*normdist1[itsindex]);
+ }
+ else{
+ if (bestConst&&bestConst->fN+bestConst->fNDeadZone>4.5){
+ dist[itsindex] = bestConst->fD[0];
+ norm[itsindex] = bestConst->fDnorm[0];
+ normdist0[itsindex] = TMath::Abs(bestConst->fD[0]/norm[itsindex]);
+ normdist1[itsindex] = TMath::Abs(bestConst->fD[0]/norm[itsindex]);
+ normdist[itsindex] = TMath::Sqrt(normdist0[itsindex]*normdist0[itsindex]+normdist1[itsindex]*normdist1[itsindex]);
+ }else{
+ dist[itsindex] = trackat0.fP0;
+ norm[itsindex] = TMath::Sqrt(trackat0.fC00);
+ normdist0[itsindex] = TMath::Abs(trackat0.fP0/norm[itsindex]);
+ normdist1[itsindex] = TMath::Abs((trackat0.fP1-primvertex[2])/TMath::Sqrt(trackat0.fC11));
+ normdist[itsindex] = TMath::Sqrt(normdist0[itsindex]*normdist0[itsindex]+normdist1[itsindex]*normdist1[itsindex]);
+ if (TMath::Abs(trackat0.fP3)>1.1){
+ if (normdist[itsindex]<10) forbidden[itsindex]=kTRUE;
+ }
+ }
+ }
+ //
+ //-----------------------------------------------------------
+ //Forbid primary track candidates -
+ //
+ //treetr->SetAlias("forbidden0","Tr0.fN<4&&Tr1.fN+Tr1.fNDeadZone>4.5");
+ //treetr->SetAlias("forbidden1","ND<3&&Tr1.fN+Tr1.fNDeadZone>5.5");
+ //treetr->SetAlias("forbidden2","ND<2&&Tr1.fClIndex[0]>0&&Tr1.fClIndex[0]>0");
+ //treetr->SetAlias("forbidden3","ND<1&&Tr1.fClIndex[0]>0");
+ //treetr->SetAlias("forbidden4","ND<4&&Tr1.fNormChi2[0]<2");
+ //treetr->SetAlias("forbidden5","ND<5&&Tr1.fNormChi2[0]<1");
+ //-----------------------------------------------------------
+ if (bestConst){
+ if (bestLong->fN<4 && bestConst->fN+bestConst->fNDeadZone>4.5) forbidden[itsindex]=kTRUE;
+ if (normdist[itsindex]<3 && bestConst->fN+bestConst->fNDeadZone>5.5) forbidden[itsindex]=kTRUE;
+ if (normdist[itsindex]<2 && bestConst->fClIndex[0]>0 && bestConst->fClIndex[1]>0 ) forbidden[itsindex]=kTRUE;
+ if (normdist[itsindex]<1 && bestConst->fClIndex[0]>0) forbidden[itsindex]=kTRUE;
+ if (normdist[itsindex]<4 && bestConst->fNormChi2[0]<2) forbidden[itsindex]=kTRUE;
+ if (normdist[itsindex]<5 && bestConst->fNormChi2[0]<1) forbidden[itsindex]=kTRUE;
+ if (bestConst->fNormChi2[0]<2.5) {
+ minPointAngle[itsindex]= 0.9999;
+ maxr[itsindex] = 10;
+ }
+ }
+ //
+ //forbid daughter kink candidates
+ //
+ if (esdtrack->GetKinkIndex(0)>0) forbidden[itsindex] = kTRUE;
+ Bool_t isElectron = kTRUE;
+ Double_t pid[5];
+ esdtrack->GetESDpid(pid);
+ for (Int_t i=1;i<5;i++){
+ if (pid[0]<pid[i]) isElectron= kFALSE;
+ }
+ if (isElectron){
+ forbidden[itsindex]=kFALSE;
+ normdist[itsindex]+=4;
+ }
+ //
+ // Causality cuts in TPC volume
+ //
+ if (esdtrack->GetTPCdensity(0,10) >0.6) maxr[itsindex] = TMath::Min(Float_t(110),maxr[itsindex]);
+ if (esdtrack->GetTPCdensity(10,30)>0.6) maxr[itsindex] = TMath::Min(Float_t(120),maxr[itsindex]);
+ if (esdtrack->GetTPCdensity(20,40)>0.6) maxr[itsindex] = TMath::Min(Float_t(130),maxr[itsindex]);
+ if (esdtrack->GetTPCdensity(30,50)>0.6) maxr[itsindex] = TMath::Min(Float_t(140),maxr[itsindex]);
+ //
+ if (esdtrack->GetTPCdensity(0,60)<0.4&&bestLong->fN<3) minr[itsindex]=100;
+ //
+ //
+ if (0){
+ cstream<<"Track"<<
+ "Tr0.="<<best<<
+ "Tr1.="<<((bestConst)? bestConst:dummy)<<
+ "Tr2.="<<bestLong<<
+ "Tr3.="<<&trackat0<<
+ "Esd.="<<esdtrack<<
+ "Dist="<<dist[itsindex]<<
+ "ND0="<<normdist0[itsindex]<<
+ "ND1="<<normdist1[itsindex]<<
+ "ND="<<normdist[itsindex]<<
+ "Pz="<<primvertex[2]<<
+ "Forbid="<<forbidden[itsindex]<<
+ "\n";
+ //
+ }
+ trackarray.AddAt(best,itsindex);
+ trackarrayc.AddAt(bestConst,itsindex);
+ trackarrayl.AddAt(bestLong,itsindex);
+ new (&helixes[itsindex]) AliHelix(*best);
+ }
//
- //find intersection linear
//
- Double_t phase[2][2],radius[2];
- Int_t points = dhelix1.GetRPHIintersections(mhelix, phase, radius);
- Double_t delta1=10000,delta2=10000;
-
- if (points>0){
- dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
- dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
- dhelix1.LinearDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
- }
- if (points==2){
- dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
- dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
- dhelix1.LinearDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
- }
- distance1 = TMath::Min(delta1,delta2);
- vertex->fDist1 = TMath::Sqrt(distance1);
//
- //find intersection parabolic
+ // first iterration of V0 finder
//
- points = dhelix1.GetRPHIintersections(mhelix, phase, radius);
- delta1=10000,delta2=10000;
-
- if (points>0){
- dhelix1.ParabolicDCA(mhelix,phase[0][0],phase[0][1],radius[0],delta1);
- }
- if (points==2){
- dhelix1.ParabolicDCA(mhelix,phase[1][0],phase[1][1],radius[1],delta2);
- }
-
- distance2 = TMath::Min(delta1,delta2);
- vertex->fDist2 = TMath::Sqrt(distance2);
- if (distance2<100){
- if (delta1<delta2){
- //get V0 info
- dhelix1.Evaluate(phase[0][0],vertex->fXr);
- dhelix1.GetMomentum(phase[0][0],vertex->fPdr);
- mhelix.GetMomentum(phase[0][1],vertex->fPm);
- dhelix1.GetAngle(phase[0][0],mhelix,phase[0][1],vertex->fAngle);
- vertex->fRr = TMath::Sqrt(radius[0]);
- }
- else{
- dhelix1.Evaluate(phase[1][0],vertex->fXr);
- dhelix1.GetMomentum(phase[1][0],vertex->fPdr);
- mhelix.GetMomentum(phase[1][1],vertex->fPm);
- dhelix1.GetAngle(phase[1][0],mhelix,phase[1][1],vertex->fAngle);
- vertex->fRr = TMath::Sqrt(radius[1]);
+ for (Int_t iesd0=0;iesd0<ntracks;iesd0++){
+ Int_t itrack0 = itsmap[iesd0];
+ if (forbidden[itrack0]) continue;
+ AliITStrackMI * btrack0 = (AliITStrackMI*)trackarray.At(itrack0);
+ if (!btrack0) continue;
+ if (btrack0->fP4>0) continue;
+ AliITStrackMI *trackc0 = (AliITStrackMI*)trackarrayc.At(itrack0);
+ //
+ for (Int_t iesd1=0;iesd1<ntracks;iesd1++){
+ Int_t itrack1 = itsmap[iesd1];
+ if (forbidden[itrack1]) continue;
+
+ AliITStrackMI * btrack1 = (AliITStrackMI*)trackarray.At(itrack1);
+ if (!btrack1) continue;
+ if (btrack1->fP4<0) continue;
+ Bool_t isGold = kFALSE;
+ if (TMath::Abs(TMath::Abs(btrack0->GetLabel())-TMath::Abs(btrack1->GetLabel()))==1){
+ isGold = kTRUE;
+ }
+ AliITStrackMI *trackc1 = (AliITStrackMI*)trackarrayc.At(itrack1);
+ AliHelix &h1 = helixes[itrack0];
+ AliHelix &h2 = helixes[itrack1];
+ //
+ // find linear distance
+ Double_t rmin =0;
+ //
+ //
+ //
+ Double_t phase[2][2],radius[2];
+ Int_t points = h1.GetRPHIintersections(h2, phase, radius);
+ if (points==0) continue;
+ Double_t delta[2]={1000,1000};
+ rmin = radius[0];
+ h1.ParabolicDCA(h2,phase[0][0],phase[0][1],radius[0],delta[0]);
+ if (points==2){
+ if (radius[1]<rmin) rmin = radius[1];
+ h1.ParabolicDCA(h2,phase[1][0],phase[1][1],radius[1],delta[1]);
+ }
+ rmin = TMath::Sqrt(rmin);
+ Double_t distance = 0;
+ Double_t radiusC = 0;
+ Int_t iphase = 0;
+ if (delta[0]<delta[1]){
+ distance = TMath::Sqrt(delta[0]);
+ radiusC = TMath::Sqrt(radius[0]);
+ }else{
+ distance = TMath::Sqrt(delta[1]);
+ radiusC = TMath::Sqrt(radius[1]);
+ iphase=1;
+ }
+ if (radiusC<TMath::Max(minr[itrack0],minr[itrack1])) continue;
+ if (radiusC>TMath::Min(maxr[itrack0],maxr[itrack1])) continue;
+ Float_t maxDist = TMath::Min(kMaxDist,Float_t(kMaxDist0+radiusC*kMaxDist1));
+ if (distance>maxDist) continue;
+ Float_t pointAngle = h1.GetPointAngle(h2,phase[iphase],primvertex);
+ if (pointAngle<TMath::Max(minPointAngle[itrack0],minPointAngle[itrack1])) continue;
+ //
+ //
+ // Double_t distance = TestV0(h1,h2,pvertex,rmin);
+ //
+ // if (distance>maxDist) continue;
+ // if (pvertex->GetRr()<kMinR) continue;
+ // if (pvertex->GetRr()>kMaxR) continue;
+ AliITStrackMI * track0=btrack0;
+ AliITStrackMI * track1=btrack1;
+ // if (pvertex->GetRr()<3.5){
+ if (radiusC<3.5){
+ //use longest tracks inside the pipe
+ track0 = (AliITStrackMI*)trackarrayl.At(itrack0);
+ track1 = (AliITStrackMI*)trackarrayl.At(itrack1);
+ }
+ //
+ //
+ pvertex->SetM(*track0);
+ pvertex->SetP(*track1);
+ pvertex->Update(primvertex);
+ if (pvertex->GetRr()<kMinR) continue;
+ if (pvertex->GetRr()>kMaxR) continue;
+ if (pvertex->GetPointAngle()<kMinPointAngle) continue;
+ if (pvertex->GetDist2()>maxDist) continue;
+ pvertex->SetLab(0,track0->GetLabel());
+ pvertex->SetLab(1,track1->GetLabel());
+
+ //
+ AliITStrackMI * htrackc0 = trackc0 ? trackc0:dummy;
+ AliITStrackMI * htrackc1 = trackc1 ? trackc1:dummy;
+
+ //
+ //
+ TObjArray * array0b = (TObjArray*)fBestHypothesys.At(itrack0);
+ if (!array0b&&pvertex->GetRr()<40 && TMath::Abs(track0->fP3)<1.1)
+ FollowProlongationTree((AliITStrackMI*)fOriginal.At(itrack0),itrack0, kFALSE);
+ TObjArray * array1b = (TObjArray*)fBestHypothesys.At(itrack1);
+ if (!array1b&&pvertex->GetRr()<40 && TMath::Abs(track1->fP3)<1.1)
+ FollowProlongationTree((AliITStrackMI*)fOriginal.At(itrack1),itrack1, kFALSE);
+ //
+ AliITStrackMI * track0b = (AliITStrackMI*)fOriginal.At(itrack0);
+ AliITStrackMI * track1b = (AliITStrackMI*)fOriginal.At(itrack1);
+ AliITStrackMI * track0l = (AliITStrackMI*)fOriginal.At(itrack0);
+ AliITStrackMI * track1l = (AliITStrackMI*)fOriginal.At(itrack1);
+
+ Float_t minchi2before0=16;
+ Float_t minchi2before1=16;
+ Float_t minchi2after0 =16;
+ Float_t minchi2after1 =16;
+ Int_t maxLayer = GetNearestLayer(pvertex->GetXrp());
+
+ if (array0b) for (Int_t i=0;i<5;i++){
+ // best track after vertex
+ AliITStrackMI * btrack = (AliITStrackMI*)array0b->At(i);
+ if (!btrack) continue;
+ if (btrack->fN>track0l->fN) track0l = btrack;
+ if (btrack->fX<pvertex->GetRr()-2) {
+ if (maxLayer>i+2 && btrack->fN==(6-i)&&i<2){
+ Float_t sumchi2= 0;
+ Float_t sumn = 0;
+ for (Int_t ilayer=i;ilayer<maxLayer;ilayer++){
+ sumn+=1.;
+ if (!btrack->fClIndex[ilayer]){
+ sumchi2+=25;
+ continue;
+ }else{
+ sumchi2+=btrack->fDy[ilayer]*btrack->fDy[ilayer]/(btrack->fSigmaY[ilayer]*btrack->fSigmaY[ilayer]);
+ sumchi2+=btrack->fDz[ilayer]*btrack->fDz[ilayer]/(btrack->fSigmaZ[ilayer]*btrack->fSigmaZ[ilayer]);
+ }
+ }
+ sumchi2/=sumn;
+ if (sumchi2<minchi2before0) minchi2before0=sumchi2;
+ }
+ continue; //safety space - Geo manager will give exact layer
+ }
+ track0b = btrack;
+ minchi2after0 = btrack->fNormChi2[i];
+ break;
+ }
+ if (array1b) for (Int_t i=0;i<5;i++){
+ // best track after vertex
+ AliITStrackMI * btrack = (AliITStrackMI*)array1b->At(i);
+ if (!btrack) continue;
+ if (btrack->fN>track1l->fN) track1l = btrack;
+ if (btrack->fX<pvertex->GetRr()-2){
+ if (maxLayer>i+2&&btrack->fN==(6-i)&&(i<2)){
+ Float_t sumchi2= 0;
+ Float_t sumn = 0;
+ for (Int_t ilayer=i;ilayer<maxLayer;ilayer++){
+ sumn+=1.;
+ if (!btrack->fClIndex[ilayer]){
+ sumchi2+=30;
+ continue;
+ }else{
+ sumchi2+=btrack->fDy[ilayer]*btrack->fDy[ilayer]/(btrack->fSigmaY[ilayer]*btrack->fSigmaY[ilayer]);
+ sumchi2+=btrack->fDz[ilayer]*btrack->fDz[ilayer]/(btrack->fSigmaZ[ilayer]*btrack->fSigmaZ[ilayer]);
+ }
+ }
+ sumchi2/=sumn;
+ if (sumchi2<minchi2before1) minchi2before1=sumchi2;
+ }
+ continue; //safety space - Geo manager will give exact layer
+ }
+ track1b = btrack;
+ minchi2after1 = btrack->fNormChi2[i];
+ break;
+ }
+ //
+ // position resolution - used for DCA cut
+ Float_t sigmad = track0b->fC00+track0b->fC11+track1b->fC00+track1b->fC11+
+ (track0b->fX-pvertex->GetRr())*(track0b->fX-pvertex->GetRr())*(track0b->fC22+track0b->fC33)+
+ (track1b->fX-pvertex->GetRr())*(track1b->fX-pvertex->GetRr())*(track1b->fC22+track1b->fC33);
+ sigmad =TMath::Sqrt(sigmad)+0.04;
+ if (pvertex->GetRr()>50){
+ Double_t cov0[15],cov1[15];
+ track0b->fESDtrack->GetInnerExternalCovariance(cov0);
+ track1b->fESDtrack->GetInnerExternalCovariance(cov1);
+ sigmad = cov0[0]+cov0[2]+cov1[0]+cov1[2]+
+ (80.-pvertex->GetRr())*(80.-pvertex->GetRr())*(cov0[5]+cov0[9])+
+ (80.-pvertex->GetRr())*(80.-pvertex->GetRr())*(cov1[5]+cov1[9]);
+ sigmad =TMath::Sqrt(sigmad)+0.05;
+ }
+ //
+ AliESDV0MI vertex2;
+ vertex2.SetM(*track0b);
+ vertex2.SetP(*track1b);
+ vertex2.Update(primvertex);
+ if (vertex2.GetDist2()<=pvertex->GetDist2()&&(vertex2.GetPointAngle()>=pvertex->GetPointAngle())){
+ pvertex->SetM(*track0b);
+ pvertex->SetP(*track1b);
+ pvertex->Update(primvertex);
+ }
+ pvertex->SetDistSigma(sigmad);
+ pvertex->SetDistNorm(pvertex->GetDist2()/sigmad);
+ //
+ // define likelihhod and causalities
+ //
+ Float_t pa0=1, pa1=1, pb0=0.26, pb1=0.26;
+ if (maxLayer<2){
+ if (pvertex->GetAnglep()[2]>0.2){
+ pb0 = TMath::Exp(-TMath::Min(normdist[itrack0],Float_t(16.))/12.);
+ pb1 = TMath::Exp(-TMath::Min(normdist[itrack1],Float_t(16.))/12.);
+ }
+ pvertex->SetChi2Before(normdist[itrack0]);
+ pvertex->SetChi2After(normdist[itrack1]);
+ pvertex->SetNAfter(0);
+ pvertex->SetNBefore(0);
+ }else{
+ pvertex->SetChi2Before(minchi2before0);
+ pvertex->SetChi2After(minchi2before1);
+ if (pvertex->GetAnglep()[2]>0.2){
+ pb0 = TMath::Exp(-TMath::Min(minchi2before0,Float_t(16))/12.);
+ pb1 = TMath::Exp(-TMath::Min(minchi2before1,float_t(16))/12.);
+ }
+ pvertex->SetNAfter(maxLayer);
+ pvertex->SetNBefore(maxLayer);
+ }
+ if (pvertex->GetRr()<90){
+ pa0 *= TMath::Min(track0->fESDtrack->GetTPCdensity(0,60),Float_t(1.));
+ pa1 *= TMath::Min(track1->fESDtrack->GetTPCdensity(0,60),Float_t(1.));
+ }
+ if (pvertex->GetRr()<20){
+ pa0 *= (0.2+TMath::Exp(-TMath::Min(minchi2after0,Float_t(16))/8.))/1.2;
+ pa1 *= (0.2+TMath::Exp(-TMath::Min(minchi2after1,Float_t(16))/8.))/1.2;
+ }
+ //
+ pvertex->SetCausality(pb0,pb1,pa0,pa1);
+ //
+ // Likelihood calculations - apply cuts
+ //
+ Bool_t v0OK = kTRUE;
+ Float_t p12 = pvertex->GetParamP()->GetParameter()[4]*pvertex->GetParamP()->GetParameter()[4];
+ p12 += pvertex->GetParamM()->GetParameter()[4]*pvertex->GetParamM()->GetParameter()[4];
+ p12 = TMath::Sqrt(p12); // "mean" momenta
+ Float_t sigmap0 = 0.0001+0.001/(0.1+pvertex->GetRr());
+ Float_t sigmap = 0.5*sigmap0*(0.6+0.4*p12); // "resolution: of point angle - as a function of radius and momenta
+
+ Float_t CausalityA = (1.0-pvertex->GetCausalityP()[0])*(1.0-pvertex->GetCausalityP()[1]);
+ Float_t CausalityB = TMath::Sqrt(TMath::Min(pvertex->GetCausalityP()[2],Float_t(0.7))*
+ TMath::Min(pvertex->GetCausalityP()[3],Float_t(0.7)));
+ //
+ Float_t Likelihood0 = (TMath::Exp(-pvertex->GetDistNorm())+0.1) *(pvertex->GetDist2()<0.5)*(pvertex->GetDistNorm()<5);
+
+ Float_t Likelihood1 = TMath::Exp(-(1.0001-pvertex->GetPointAngle())/sigmap)+
+ 0.4*TMath::Exp(-(1.0001-pvertex->GetPointAngle())/(4.*sigmap))+
+ 0.4*TMath::Exp(-(1.0001-pvertex->GetPointAngle())/(8.*sigmap))+
+ 0.1*TMath::Exp(-(1.0001-pvertex->GetPointAngle())/0.01);
+ //
+ if (CausalityA<kCausality0Cut) v0OK = kFALSE;
+ if (TMath::Sqrt(Likelihood0*Likelihood1)<kLikelihood01Cut) v0OK = kFALSE;
+ if (Likelihood1<kLikelihood1Cut) v0OK = kFALSE;
+ if (TMath::Power(Likelihood0*Likelihood1*CausalityB,0.33)<kCombinedCut) v0OK = kFALSE;
+
+ //
+ //
+ if (kTRUE){
+ Bool_t gold = TMath::Abs(TMath::Abs(track0->GetLabel())-TMath::Abs(track1->GetLabel()))==1;
+ cstream<<"It0"<<
+ "Tr0.="<<track0<< //best without constrain
+ "Tr1.="<<track1<< //best without constrain
+ "Tr0B.="<<track0b<< //best without constrain after vertex
+ "Tr1B.="<<track1b<< //best without constrain after vertex
+ "Tr0C.="<<htrackc0<< //best with constrain if exist
+ "Tr1C.="<<htrackc1<< //best with constrain if exist
+ "Tr0L.="<<track0l<< //longest best
+ "Tr1L.="<<track1l<< //longest best
+ "Esd0.="<<track0->fESDtrack<< // esd track0 params
+ "Esd1.="<<track1->fESDtrack<< // esd track1 params
+ "V0.="<<pvertex<< //vertex properties
+ "V0b.="<<&vertex2<< //vertex properties at "best" track
+ "ND0="<<normdist[itrack0]<< //normalize distance for track0
+ "ND1="<<normdist[itrack1]<< //normalize distance for track1
+ "Gold="<<gold<< //
+ // "RejectBase="<<rejectBase<< //rejection in First itteration
+ "OK="<<v0OK<<
+ "rmin="<<rmin<<
+ "sigmad="<<sigmad<<
+ "\n";
+ }
+ //if (rejectBase) continue;
+ //
+ pvertex->SetStatus(0);
+ // if (rejectBase) {
+ // pvertex->SetStatus(-100);
+ //}
+ if (pvertex->GetPointAngle()>kMinPointAngle2) {
+ if (v0OK){
+ AliV0vertex vertexjuri(*track0,*track1);
+ vertexjuri.SetESDindexes(track0->fESDtrack->GetID(),track1->fESDtrack->GetID());
+ pvertex->SetESDindexes(track0->fESDtrack->GetID(),track1->fESDtrack->GetID());
+ event->AddV0(&vertexjuri);
+ pvertex->SetStatus(100);
+ }
+ event->AddV0MI(pvertex);
+ }
}
}
- //
//
- return vertex->fDist2;
+ //
+ // delete temporary arrays
+ //
+ delete[] minPointAngle;
+ delete[] maxr;
+ delete[] minr;
+ delete[] norm;
+ delete[] normdist;
+ delete[] normdist1;
+ delete[] normdist0;
+ delete[] dist;
+ delete[] itsmap;
+ delete[] helixes;
+ delete pvertex;
}
+
+
+
+
+
+
+
+
+
+
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff