--- /dev/null
+/****************************************************************************
+ * This macro is supposed to do reconstruction in the barrel ALICE trackers *
+ * (Make sure you have TPC digits and ITS hits before using this macro !!!) *
+ * It does the following steps (April 12, 2001): *
+ * 1) TPC cluster finding *
+ * 2) TPC track finding *
+ * 3) ITS cluster finding V2 (via fast points !) *
+ * 4) ITS track finding V2 *
+ * 5) ITS back track propagation V2 *
+ * 6) TPC back track propagation *
+ * (Origin: I.Belikov, CERN, Jouri.Belikov@cern.ch) *
+ ****************************************************************************/
+
+#ifndef __CINT__
+ #include "alles.h"
+ #include "AliMagF.h"
+ #include "AliTPCtracker.h"
+
+ #include "AliITS.h"
+ #include "AliITSgeom.h"
+ #include "AliITSRecPoint.h"
+ #include "AliITSclusterV2.h"
+ #include "AliITSsimulationFastPoints.h"
+ #include "AliITStrackerV2.h"
+#endif
+
+Int_t TPCFindClusters(const Char_t *inname, const Char_t *outname);
+Int_t TPCFindTracks(const Char_t *inname, const Char_t *outname);
+Int_t TPCSortTracks(const Char_t *inname, const Char_t *outname);
+Int_t TPCPropagateBack(const Char_t *inname, const Char_t *outname);
+
+Int_t ITSFindClusters(const Char_t *inname, const Char_t *outname);
+Int_t ITSFindTracks(const Char_t *inname, const Char_t *outname);
+Int_t ITSPropagateBack(const Char_t *inname, const Char_t *outname);
+
+
+Int_t AliBarrelReconstruction() {
+ const Char_t *TPCdigName="galice.root";
+ const Char_t *TPCclsName="AliTPCclusters.root";
+ const Char_t *TPCtrkName="AliTPCtracks.root";
+
+ const Char_t *ITSdigName="galice.root";
+ const Char_t *ITSclsName="AliITSclustersV2.root";
+ const Char_t *ITStrkName="AliITStracksV2.root";
+
+ AliKalmanTrack::SetConvConst(100/0.299792458/0.2/gAlice->Field()->Factor());
+
+// ********** Find TPC clusters *********** //
+ if (TPCFindClusters(TPCdigName,TPCclsName)) {
+ cerr<<"Failed to get TPC clusters !\n";
+ return 1;
+ }
+
+// ********** Find TPC tracks *********** //
+ if (TPCFindTracks(TPCclsName,TPCtrkName)) {
+ cerr<<"Failed to get TPC tracks !\n";
+ return 1;
+ }
+
+// ********** Sort and label TPC tracks *********** //
+ if (TPCSortTracks(TPCclsName,TPCtrkName)) {
+ cerr<<"Failed to sort TPC tracks !\n";
+ return 1;
+ }
+
+// ********** Find ITS clusters *********** //
+ if (ITSFindClusters(ITSdigName,ITSclsName)) {
+ cerr<<"Failed to get ITS clusters !\n";
+ return 1;
+ }
+
+// ********** Find ITS tracks *********** //
+ {TFile *clsFile=TFile::Open(ITSclsName);
+ if (ITSFindTracks(TPCtrkName,ITStrkName)) {
+ cerr<<"Failed to get ITS tracks !\n";
+ return 1;
+ }
+ clsFile->Close();}
+
+// ********** Back propagation of the ITS tracks *********** //
+ {TFile *clsFile=TFile::Open(ITSclsName);
+ if (ITSPropagateBack(ITStrkName,TPCtrkName)) {
+ cerr<<"Failed to propagate back the ITS tracks !\n";
+ return 1;
+ }
+ clsFile->Close();}
+
+
+// ********** Back propagation of the TPC tracks *********** //
+ {TFile *clsFile=TFile::Open(TPCclsName);
+ if (TPCPropagateBack(TPCtrkName,TPCtrkName)) {
+ cerr<<"Failed to propagate back the TPC tracks !\n";
+ return 1;
+ }
+ clsFile->Close();}
+
+ return 0;
+}
+
+
+Int_t TPCFindClusters(const Char_t *inname, const Char_t *outname) {
+ Int_t rc=0;
+ const Char_t *name="TPCFindClusters";
+ cerr<<'\n'<<name<<"...\n";
+ gBenchmark->Start(name);
+
+ TFile *out=TFile::Open(outname,"recreate");
+ TFile *in =TFile::Open(inname);
+
+ AliTPCParam *param=(AliTPCParam *)in->Get("75x40_100x60");
+ if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;}
+ AliTPCv2 tpc;
+ tpc.SetParam(param);
+
+ tpc.Digits2Clusters(out);
+
+ in->Close();
+ out->Close();
+ gBenchmark->Stop(name);
+ gBenchmark->Show(name);
+
+ return rc;
+}
+
+Int_t TPCFindTracks(const Char_t *inname, const Char_t *outname) {
+ Int_t rc=0;
+ const Char_t *name="TPCFindTracks";
+ cerr<<'\n'<<name<<"...\n";
+ gBenchmark->Start(name);
+ TFile *out=TFile::Open(outname,"recreate");
+ TFile *in =TFile::Open(inname);
+ AliTPCParam *param=(AliTPCParam *)in->Get("75x40_100x60");
+ if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;}
+
+ AliTPCtracker *tracker=new AliTPCtracker(param);
+ rc=tracker->Clusters2Tracks(0,out);
+ delete tracker;
+
+ in->Close();
+ out->Close();
+ gBenchmark->Stop(name);
+ gBenchmark->Show(name);
+
+ return rc;
+}
+
+Int_t TPCSortTracks(const Char_t *inname, const Char_t *outname) {
+ Int_t rc=0;
+ const Char_t *name="TPCSortTracks";
+ cerr<<'\n'<<name<<"...\n";
+ gBenchmark->Start(name);
+
+ TFile *out=TFile::Open(outname);
+ TFile *in =TFile::Open(inname);
+ AliTPCParam *param=(AliTPCParam *)in->Get("75x40_100x60");
+ if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;}
+
+ TObjArray tarray(10000);
+ AliTPCtrack *iotrack=0;
+ Int_t i;
+
+ out->cd();
+ TTree *tracktree=(TTree*)out->Get("TPCf");
+ TBranch *tbranch=tracktree->GetBranch("tracks");
+ Int_t nentr=(Int_t)tracktree->GetEntries();
+ for (i=0; i<nentr; i++) {
+ iotrack=new AliTPCtrack;
+ tbranch->SetAddress(&iotrack);
+ tracktree->GetEvent(i);
+ tarray.AddLast(iotrack);
+ }
+ tarray.Sort();
+ out->Close();
+
+ //assign thacks GEANT labels
+ in->cd();
+ AliTPCtracker *tracker = new AliTPCtracker(param);
+ tracker->LoadInnerSectors();
+ tracker->LoadOuterSectors();
+ for (i=0; i<nentr; i++) {
+ iotrack=(AliTPCtrack*)tarray.UncheckedAt(i);
+ tracker->CookLabel(iotrack,0.1);
+ }
+ delete tracker;
+ in->Close();
+ //end of GEANT label assignment
+
+ out=TFile::Open(outname,"recreate");
+ tracktree=new TTree("TPCf","Tree with TPC tracks");
+ tracktree->Branch("tracks","AliTPCtrack",&iotrack,32000,0);
+ for (i=0; i<nentr; i++) {
+ iotrack=(AliTPCtrack*)tarray.UncheckedAt(i);
+ tracktree->Fill();
+ }
+ tracktree->Write();
+ out->Close();
+
+ gBenchmark->Stop(name);
+ gBenchmark->Show(name);
+
+ return rc;
+}
+
+Int_t ITSFindClusters(const Char_t *inname, const Char_t *outname) {
+ Int_t rc=0;
+ const Char_t *name="ITSFindClusters";
+ cerr<<'\n'<<name<<"...\n";
+ gBenchmark->Start(name);
+ TFile *out=TFile::Open(outname,"recreate");
+ TFile *in =TFile::Open(inname,"update");
+
+ if (!(gAlice=(AliRun*)in->Get("gAlice"))) {
+ cerr<<"Can't get gAlice !\n";
+ return 1;
+ }
+ Int_t ev=0;
+ gAlice->GetEvent(ev);
+ AliITS *ITS = (AliITS*)gAlice->GetModule("ITS");
+ if (!ITS) { cerr<<"Can't get the ITS !\n"; return 1;}
+
+ gAlice->MakeTree("R"); ITS->MakeBranch("R",0);
+ //////////////// Taken from ITSHitsToFastPoints.C ///////////////////////
+ AliITSsimulationFastPoints *sim = new AliITSsimulationFastPoints();
+ for (Int_t i=0;i<3;i++) { ITS->SetSimulationModel(i,sim); }
+ Int_t nsignal=25;
+ Int_t size=-1;
+ Int_t bgr_ev=Int_t(ev/nsignal);
+ ITS->HitsToFastRecPoints(ev,bgr_ev,size," ","All"," ");
+ //////////////////////////////////////////////////////////////////////////
+
+ gAlice->GetEvent(ev);
+
+ AliITSgeom *geom=ITS->GetITSgeom();
+ out->cd();
+ geom->Write();
+
+ TClonesArray *clusters=new TClonesArray("AliITSclusterV2",10000);
+ TTree *cTree=new TTree("cTree","ITS clusters");
+ cTree->Branch("Clusters",&clusters);
+
+ TTree *pTree=gAlice->TreeR();
+ if (!pTree) { cerr<<"Can't get TreeR !\n"; return 1; }
+ TBranch *branch=pTree->GetBranch("ITSRecPoints");
+ if (!branch) { cerr<<"Can't get ITSRecPoints branch !\n"; return 1;}
+ TClonesArray *points=new TClonesArray("AliITSRecPoint",10000);
+ branch->SetAddress(&points);
+
+ TClonesArray &cl=*clusters;
+ Int_t nclusters=0;
+ Int_t nentr=(Int_t)pTree->GetEntries();
+ for (Int_t i=0; i<nentr; i++) {
+ if (!pTree->GetEvent(i)) {cTree->Fill(); continue;}
+ Int_t lay,lad,det; geom->GetModuleId(i,lay,lad,det);
+ Float_t x,y,zshift; geom->GetTrans(lay,lad,det,x,y,zshift);
+ Double_t rot[9]; geom->GetRotMatrix(lay,lad,det,rot);
+ Double_t yshift = x*rot[0] + y*rot[1];
+ Int_t ndet=(lad-1)*geom->GetNdetectors(lay) + (det-1);
+ Int_t ncl=points->GetEntriesFast();
+ nclusters+=ncl;
+ for (Int_t j=0; j<ncl; j++) {
+ AliITSRecPoint *p=(AliITSRecPoint*)points->UncheckedAt(j);
+ Float_t lp[5];
+ lp[0]=-p->GetX()-yshift; if (lay==1) lp[0]=-lp[0];
+ lp[1]=p->GetZ()+zshift;
+ lp[2]=p->GetSigmaX2();
+ lp[3]=p->GetSigmaZ2();
+ lp[4]=p->GetQ();
+ Int_t lab[6];
+ lab[0]=p->GetLabel(0);lab[1]=p->GetLabel(1);lab[2]=p->GetLabel(2);
+ lab[3]=ndet;
+
+ Int_t label=lab[0];
+ TParticle *part=(TParticle*)gAlice->Particle(label);
+ label=-3;
+ while (part->P() < 0.005) {
+ Int_t m=part->GetFirstMother();
+ if (m<0) {cerr<<"Primary momentum: "<<part->P()<<endl; break;}
+ label=m;
+ part=(TParticle*)gAlice->Particle(label);
+ }
+ if (lab[1]<0) lab[1]=label;
+ else if (lab[2]<0) lab[2]=label;
+ else cerr<<"No empty labels !\n";
+
+ new(cl[j]) AliITSclusterV2(lab,lp);
+ }
+ cTree->Fill(); clusters->Delete();
+ points->Delete();
+ }
+ cTree->Write();
+ cerr<<"Number of clusters: "<<nclusters<<endl;
+
+ delete cTree; delete clusters; delete points;
+
+ delete gAlice; gAlice=0;
+ in->Close();
+ out->Close();
+ gBenchmark->Stop(name);
+ gBenchmark->Show(name);
+
+ return rc;
+}
+
+Int_t ITSFindTracks(const Char_t *inname, const Char_t *outname) {
+ Int_t rc=0;
+ const Char_t *name="ITSFindTracks";
+ cerr<<'\n'<<name<<"...\n";
+ gBenchmark->Start(name);
+
+ AliITSgeom *geom=(AliITSgeom*)gFile->Get("AliITSgeom");
+ if (!geom) { cerr<<"can't get ITS geometry !\n"; return 1;}
+ AliITStrackerV2 tracker(geom);
+
+ TFile *out=TFile::Open(outname,"recreate");
+ TFile *in =TFile::Open(inname);
+ rc=tracker.Clusters2Tracks(in,out);
+ in->Close();
+ out->Close();
+
+ gBenchmark->Stop(name);
+ gBenchmark->Show(name);
+
+ return rc;
+}
+
+Int_t ITSPropagateBack(const Char_t *inname, const Char_t *outname) {
+ Int_t rc=0;
+ const Char_t *name="ITSPropagateBack";
+ cerr<<'\n'<<name<<"...\n";
+ gBenchmark->Start(name);
+
+ AliITSgeom *geom=(AliITSgeom*)gFile->Get("AliITSgeom");
+ if (!geom) { cerr<<"can't get ITS geometry !\n"; return 1;}
+ AliITStrackerV2 tracker(geom);
+
+ TFile *out=TFile::Open(outname,"update");
+ TFile *in =TFile::Open(inname);
+ rc=tracker.PropagateBack(in,out);
+ in->Close();
+ out->Close();
+
+ gBenchmark->Stop(name);
+ gBenchmark->Show(name);
+
+ return rc;
+}
+
+Int_t TPCPropagateBack(const Char_t *inname, const Char_t *outname) {
+ Int_t rc=0;
+ const Char_t *name="TPCPropagateBack";
+ cerr<<'\n'<<name<<"...\n";
+ gBenchmark->Start(name);
+
+ AliTPCParam *param=(AliTPCParam *)gFile->Get("75x40_100x60");
+ if (!param) {cerr<<"TPC parameters have not been found !\n"; return 1;}
+ AliTPCtracker *tracker=new AliTPCtracker(param);
+
+ TFile *out=TFile::Open(outname,"update");
+ TFile *in =TFile::Open(inname);
+ rc=tracker->PropagateBack(in,out);
+ delete tracker;
+ in->Close();
+ out->Close();
+
+ gBenchmark->Stop(name);
+ gBenchmark->Show(name);
+
+ return rc;
+}
+
/*
$Log$
+Revision 1.33 2001/04/03 12:40:43 kowal2
+Removed printouts
+
Revision 1.32 2001/03/12 17:47:36 hristov
Changes needed on Sun with CC 5.0
//-----------------------------------------------------------------
// This is a track finder.
//-----------------------------------------------------------------
- AliTPCtracker::Clusters2Tracks(fTPCParam,of);
+ AliTPCtracker tracker(fTPCParam);
+ tracker.Clusters2Tracks(gFile,of);
}
//_____________________________________________________________________________
#ifndef __CINT__
#include "alles.h"
+ #include "AliTPCtracker.h"
#endif
struct GoodTrack {
if (!cf->IsOpen()) {cerr<<"Can't open AliTPCclusters.root !\n"; return 1;}
AliTPCParam *digp= (AliTPCParam*)cf->Get("75x40_100x60");
if (!digp) { cerr<<"TPC parameters have not been found !\n"; return 2; }
+ AliTPCtracker *tracker = new AliTPCtracker(digp);
// Load clusters
- AliTPCClustersArray *ca=new AliTPCClustersArray;
- ca->Setup(digp);
- ca->SetClusterType("AliTPCcluster");
- ca->ConnectTree("Segment Tree");
- Int_t nentr=Int_t(ca->GetTree()->GetEntries());
- for (i=0; i<nentr; i++) ca->LoadEntry(i);
+ tracker->LoadInnerSectors();
+ tracker->LoadOuterSectors();
// Load tracks
TFile *tf=TFile::Open("AliTPCtracks.root");
if (!tf->IsOpen()) {cerr<<"Can't open AliTPCtracks.root !\n"; return 3;}
TObjArray tarray(2000);
- TTree *tracktree=(TTree*)tf->Get("TreeT");
+ TTree *tracktree=(TTree*)tf->Get("TPCf");
+ if (!tracktree) {cerr<<"Can't get a tree with TPC tracks !\n"; return 4;}
TBranch *tbranch=tracktree->GetBranch("tracks");
- nentr=(Int_t)tracktree->GetEntries();
+ Int_t nentr=(Int_t)tracktree->GetEntries();
AliTPCtrack *iotrack=0;
for (i=0; i<nentr; i++) {
iotrack=new AliTPCtrack;
tbranch->SetAddress(&iotrack);
tracktree->GetEvent(i);
- iotrack->CookLabel(ca);
+ tracker->CookLabel(iotrack,0.1);
tarray.AddLast(iotrack);
}
- tf->Close();
- delete ca;
+ delete tracker;
+
+ tf->Close();
cf->Close();
/////////////////////////////////////////////////////////////////////////
cerr<<"Preparing tracks for matching with the ITS...\n";
tarray.Sort();
tf=TFile::Open("AliTPCtracks.root","recreate");
- tracktree=new TTree("TreeT","Tree with TPC tracks");
+ tracktree=new TTree("TPCf","Tree with TPC tracks");
tracktree->Branch("tracks","AliTPCtrack",&iotrack,32000,0);
for (i=0; i<nentr; i++) {
iotrack=(AliTPCtrack*)tarray.UncheckedAt(i);
}
tracktree->Write();
tf->Close();
-
}
cerr<<"Number of good tracks : "<<ngood<<endl;
TH1F *hmpt=new TH1F("hmpt","Relative Pt resolution (pt>4GeV/c)",30,-60,60);
hmpt->SetFillColor(6);
- TH1F *hgood=new TH1F("hgood","Good tracks",30,0.1,6.1);
- TH1F *hfound=new TH1F("hfound","Found tracks",30,0.1,6.1);
- TH1F *hfake=new TH1F("hfake","Fake tracks",30,0.1,6.1);
- TH1F *hg=new TH1F("hg","",30,0.1,6.1); //efficiency for good tracks
+ AliTPCtrack *trk=(AliTPCtrack*)tarray.UncheckedAt(0);
+ Double_t pmin=0.1*(100/0.299792458/0.2/trk->GetConvConst());
+ Double_t pmax=6.0+pmin;
+
+ TH1F *hgood=new TH1F("hgood","Good tracks",30,pmin,pmax);
+ TH1F *hfound=new TH1F("hfound","Found tracks",30,pmin,pmax);
+ TH1F *hfake=new TH1F("hfake","Fake tracks",30,pmin,pmax);
+ TH1F *hg=new TH1F("hg","",30,pmin,pmax); //efficiency for good tracks
hg->SetLineColor(4); hg->SetLineWidth(2);
- TH1F *hf=new TH1F("hf","Efficiency for fake tracks",30,0.1,6.1);
+ TH1F *hf=new TH1F("hf","Efficiency for fake tracks",30,pmin,pmax);
hf->SetFillColor(1); hf->SetFillStyle(3013); hf->SetLineWidth(2);
TH1F *he =new TH1F("he","dE/dX for pions with 0.4<p<0.5 GeV/c",50,0.,100.);
Float_t ptg=
TMath::Sqrt(gt[ngood].px*gt[ngood].px + gt[ngood].py*gt[ngood].py);
+ if (ptg<pmin) continue;
+
hgood->Fill(ptg);
AliTPCtrack *track=0;
hg->SetXTitle("Pt (GeV/c)");
hg->Draw();
- TLine *line1 = new TLine(0.1,1.0,6.1,1.0); line1->SetLineStyle(4);
+ TLine *line1 = new TLine(pmin,1.0,pmax,1.0); line1->SetLineStyle(4);
line1->Draw("same");
- TLine *line2 = new TLine(0.1,0.9,6.1,0.9); line2->SetLineStyle(4);
+ TLine *line2 = new TLine(pmin,0.9,pmax,0.9); line2->SetLineStyle(4);
line2->Draw("same");
hf->SetFillColor(1);
Int_t good_tracks(GoodTrack *gt, Int_t max) {
Int_t nt=0;
- //TFile *file=TFile::Open("rfio:galice.root");
TFile *file=TFile::Open("galice.root");
if (!file->IsOpen()) {cerr<<"Can't open galice.root !\n"; exit(4);}
+#ifndef __CINT__
+ #include <iostream.h>
+ #include "AliRun.h"
+ #include "AliTPCv1.h"
+ #include "AliTPCv2.h"
+ #include "AliTPCParam.h"
+
+ #include "TFile.h"
+ #include "TStopwatch.h"
+#endif
+
Int_t AliTPCFindClusters() {
TFile *out=TFile::Open("AliTPCclusters.root","new");
if (!out->IsOpen()) {cerr<<"Delete old AliTPCclusters.root !\n"; return 1;}
+#ifndef __CINT__
+ #include <iostream.h>
+ #include "AliTPCtracker.h"
+
+ #include "TFile.h"
+ #include "TStopwatch.h"
+#endif
+
Int_t AliTPCFindTracks() {
cerr<<"Looking for tracks...\n";
TFile *in=TFile::Open("AliTPCclusters.root");
if (!in->IsOpen()) {cerr<<"Can't open AliTPCclusters.root !\n"; return 2;}
- AliTPCv2 TPC;
- AliTPCParam *digp= (AliTPCParam*)in->Get("75x40_100x60");
- if (!digp) {cerr<<"TPC parameters have not been found !\n"; return 3;}
- TPC.SetParam(digp);
-
+ AliTPCParam *par=(AliTPCParam*)in->Get("75x40_100x60");
+ if (!par) {cerr<<"Can't get TPC parameters !\n"; return 3;}
+
TStopwatch timer;
- TPC.Clusters2Tracks(out);
+ AliTPCtracker *tracker = new AliTPCtracker(par);
+ Int_t rc=tracker->Clusters2Tracks(0,out);
+ delete tracker;
timer.Stop(); timer.Print();
-
+
in->Close();
out->Close();
- return 0;
+ return rc;
}
--- /dev/null
+#ifndef __CINT__
+ #include <iostream.h>
+ #include "AliTPCtracker.h"
+
+ #include "TFile.h"
+ #include "TStopwatch.h"
+#endif
+
+Int_t AliTPCPropagateBack() {
+ cerr<<"Propagating tracks back through the TPC...\n";
+
+ TFile *in=TFile::Open("AliTPCtracks.root");
+ if (!in->IsOpen()) {
+ cerr<<"Can't open AliTPCtracks.root !\n";
+ return 1;
+ }
+ TFile *out=TFile::Open("AliTPCBackTracks.root","new");
+ if (!out->IsOpen()) {
+ cerr<<"Delete old AliTPCBackTracks.root !\n"; return 2;
+ }
+ TFile *file=TFile::Open("AliTPCclusters.root");
+ if (!file->IsOpen()) {
+ cerr<<"Can't open AliTPCclusters.root !\n";return 3;
+ }
+ AliTPCParam *param=(AliTPCParam*)file->Get("75x40_100x60");
+
+ TStopwatch timer;
+ AliTPCtracker *tracker=new AliTPCtracker(param);
+ Int_t rc=tracker->PropagateBack(in,out);
+ delete tracker;
+ timer.Stop(); timer.Print();
+
+ file->Close();
+ in->Close();
+ out->Close();
+
+ return rc;
+}
--- /dev/null
+#ifndef ALITPCRECO_H
+#define ALITPCRECO_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//-------------------------------------------------------------------------
+// TPC reconstruction name space
+//
+// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
+//-------------------------------------------------------------------------
+#include <Rtypes.h>
+
+//namespace AliTPCreco {
+ const Int_t kMaxClusterPerRow=2500;
+ const Int_t kRowsToSkip=10;
+ const Int_t kMaxRow=200;
+
+ const Double_t kMaxCHI2=12.;
+ const Double_t kMaxROAD=30.;
+//}
+
+//using namespace AliTPCreco;
+
+#endif
#include <iostream.h>
-#include "AliCluster.h"
#include "AliTPCtrack.h"
#include "AliTPCcluster.h"
#include "AliTPCClustersRow.h"
#include "AliTPCClustersArray.h"
ClassImp(AliTPCtrack)
+
//_________________________________________________________________________
AliTPCtrack::AliTPCtrack(UInt_t index, const Double_t xx[5],
const Double_t cc[15], Double_t xref, Double_t alpha) : AliKalmanTrack() {
fIndex[0]=index;
SetNumberOfClusters(1);
+
+}
+
+//_____________________________________________________________________________
+AliTPCtrack::AliTPCtrack(const AliKalmanTrack& t,Double_t alpha) {
+ //-----------------------------------------------------------------
+ // Conversion AliKalmanTrack -> AliTPCtrack.
+ //-----------------------------------------------------------------
+ SetLabel(t.GetLabel());
+ SetChi2(0.);
+ SetNumberOfClusters(0);
+ //SetConvConst(t.GetConvConst());
+
+ fdEdx = 0.;
+ fAlpha = alpha;
+ if (fAlpha < -TMath::Pi()) fAlpha += 2*TMath::Pi();
+ else if (fAlpha >= TMath::Pi()) fAlpha -= 2*TMath::Pi();
+
+ //Conversion of the track parameters
+ Double_t x,p[5]; t.GetExternalParameters(x,p);
+ fX=x; x=GetConvConst();
+ fP0=p[0];
+ fP1=p[1];
+ fP3=p[3];
+ fP4=p[4]/x;
+ fP2=fP4*fX - p[2];
+
+ //Conversion of the covariance matrix
+ Double_t c[15]; t.GetExternalCovariance(c);
+ c[10]/=x; c[11]/=x; c[12]/=x; c[13]/=x; c[14]/=x*x;
+
+ Double_t c22=fX*fX*c[14] - 2*fX*c[12] + c[5];
+ Double_t c32=fX*c[13] - c[8];
+ Double_t c20=fX*c[10] - c[3], c21=fX*c[11] - c[4], c42=fX*c[14] - c[12];
+
+ fC00=c[0 ];
+ fC10=c[1 ]; fC11=c[2 ];
+ fC20=c20; fC21=c21; fC22=c22;
+ fC30=c[6 ]; fC31=c[7 ]; fC32=c32; fC33=c[9 ];
+ fC40=c[10]; fC41=c[11]; fC42=c42; fC43=c[13]; fC44=c[14];
+
}
//_____________________________________________________________________________
//-------------------------------------------------------------------------
Double_t a=GetConvConst();
- Double_t c22=fX*fX*fC33-2*fX*fC32+fC22;
- Double_t c42=fX*fC43-fC42;
- Double_t c20=fX*fC30-fC20, c21=fX*fC31-fC21, c32=fX*fC33-fC32;
+ Double_t c22=fX*fX*fC44-2*fX*fC42+fC22;
+ Double_t c32=fX*fC43-fC32;
+ Double_t c20=fX*fC40-fC20, c21=fX*fC41-fC21, c42=fX*fC44-fC42;
cc[0 ]=fC00;
cc[1 ]=fC10; cc[2 ]=fC11;
cc[3 ]=c20; cc[4 ]=c21; cc[5 ]=c22;
- cc[6 ]=fC40; cc[7 ]=fC41; cc[8 ]=c42; cc[9 ]=fC44;
- cc[10]=fC30*a; cc[11]=fC31*a; cc[12]=c32*a; cc[13]=fC43*a; cc[14]=fC33*a*a;
+ cc[6 ]=fC30; cc[7 ]=fC31; cc[8 ]=c32; cc[9 ]=fC33;
+ cc[10]=fC40*a; cc[11]=fC41*a; cc[12]=c42*a; cc[13]=fC43*a; cc[14]=fC44*a*a;
}
//-----------------------------------------------------------------
// This function propagates a track to a reference plane x=xk.
//-----------------------------------------------------------------
- if (TMath::Abs(fP3*xk - fP2) >= 0.99999) {
+ if (TMath::Abs(fP4*xk - fP2) >= 0.99999) {
Int_t n=GetNumberOfClusters();
if (n>4) cerr<<n<<" AliTPCtrack warning: Propagation failed !\n";
return 0;
}
Double_t x1=fX, x2=x1+(xk-x1), dx=x2-x1, y1=fP0, z1=fP1;
- Double_t c1=fP3*x1 - fP2, r1=sqrt(1.- c1*c1);
- Double_t c2=fP3*x2 - fP2, r2=sqrt(1.- c2*c2);
+ Double_t c1=fP4*x1 - fP2, r1=sqrt(1.- c1*c1);
+ Double_t c2=fP4*x2 - fP2, r2=sqrt(1.- c2*c2);
fP0 += dx*(c1+c2)/(r1+r2);
- fP1 += dx*(c1+c2)/(c1*r2 + c2*r1)*fP4;
+ fP1 += dx*(c1+c2)/(c1*r2 + c2*r1)*fP3;
//f = F - 1
Double_t rr=r1+r2, cc=c1+c2, xx=x1+x2;
Double_t f02=-dx*(2*rr + cc*(c1/r1 + c2/r2))/(rr*rr);
- Double_t f03= dx*(rr*xx + cc*(c1*x1/r1+c2*x2/r2))/(rr*rr);
+ Double_t f04= dx*(rr*xx + cc*(c1*x1/r1+c2*x2/r2))/(rr*rr);
Double_t cr=c1*r2+c2*r1;
- Double_t f12=-dx*fP4*(2*cr + cc*(c2*c1/r1-r1 + c1*c2/r2-r2))/(cr*cr);
- Double_t f13=dx*fP4*(cr*xx-cc*(r1*x2-c2*c1*x1/r1+r2*x1-c1*c2*x2/r2))/(cr*cr);
- Double_t f14= dx*cc/cr;
+ Double_t f12=-dx*fP3*(2*cr + cc*(c2*c1/r1-r1 + c1*c2/r2-r2))/(cr*cr);
+ Double_t f13= dx*cc/cr;
+ Double_t f14=dx*fP3*(cr*xx-cc*(r1*x2-c2*c1*x1/r1+r2*x1-c1*c2*x2/r2))/(cr*cr);
//b = C*ft
- Double_t b00=f02*fC20 + f03*fC30, b01=f12*fC20 + f13*fC30 + f14*fC40;
- Double_t b10=f02*fC21 + f03*fC31, b11=f12*fC21 + f13*fC31 + f14*fC41;
- Double_t b20=f02*fC22 + f03*fC32, b21=f12*fC22 + f13*fC32 + f14*fC42;
- Double_t b30=f02*fC32 + f03*fC33, b31=f12*fC32 + f13*fC33 + f14*fC43;
- Double_t b40=f02*fC42 + f03*fC43, b41=f12*fC42 + f13*fC43 + f14*fC44;
+ Double_t b00=f02*fC20 + f04*fC40, b01=f12*fC20 + f14*fC40 + f13*fC30;
+ Double_t b10=f02*fC21 + f04*fC41, b11=f12*fC21 + f14*fC41 + f13*fC31;
+ Double_t b20=f02*fC22 + f04*fC42, b21=f12*fC22 + f14*fC42 + f13*fC32;
+ Double_t b30=f02*fC32 + f04*fC43, b31=f12*fC32 + f14*fC43 + f13*fC33;
+ Double_t b40=f02*fC42 + f04*fC44, b41=f12*fC42 + f14*fC44 + f13*fC43;
//a = f*b = f*C*ft
- Double_t a00=f02*b20+f03*b30,a01=f02*b21+f03*b31,a11=f12*b21+f13*b31+f14*b41;
+ Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a11=f12*b21+f14*b41+f13*b31;
//F*C*Ft = C + (a + b + bt)
fC00 += a00 + 2*b00;
Double_t theta2=14.1*14.1/(beta2*p2*1e6)*d/x0*rho;
//Double_t theta2=1.0259e-6*10*10/20/(beta2*p2)*d*rho;
- Double_t ey=fP3*fX - fP2, ez=fP4;
- Double_t xz=fP3*ez, zz1=ez*ez+1, xy=fP2+ey;
+ Double_t ey=fP4*fX - fP2, ez=fP3;
+ Double_t xz=fP4*ez, zz1=ez*ez+1, xy=fP2+ey;
- fC33 += xz*xz*theta2;
- fC32 += xz*ez*xy*theta2;
- fC43 += xz*zz1*theta2;
fC22 += (2*ey*ez*ez*fP2+1-ey*ey+ez*ez+fP2*fP2*ez*ez)*theta2;
- fC42 += ez*zz1*xy*theta2;
- fC44 += zz1*zz1*theta2;
+ fC32 += ez*zz1*xy*theta2;
+ fC33 += zz1*zz1*theta2;
+ fC42 += xz*ez*xy*theta2;
+ fC43 += xz*zz1*theta2;
+ fC44 += xz*xz*theta2;
//Energy losses************************
Double_t dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2)) - beta2)*d*rho;
if (x1 < x2) dE=-dE;
- cc=fP3;
- fP3*=(1.- sqrt(p2+pm*pm)/p2*dE);
- fP2+=fX*(fP3-cc);
+ cc=fP4;
+ fP4*=(1.- sqrt(p2+pm*pm)/p2*dE);
+ fP2+=fX*(fP4-cc);
return 1;
}
//-----------------------------------------------------------------
// This function propagates tracks to the "vertex".
//-----------------------------------------------------------------
- Double_t c=fP3*fX - fP2;
- Double_t tgf=-fP2/(fP3*fP0 + sqrt(1-c*c));
+ Double_t c=fP4*fX - fP2;
+ Double_t tgf=-fP2/(fP4*fP0 + sqrt(1-c*c));
Double_t snf=tgf/sqrt(1.+ tgf*tgf);
- Double_t xv=(fP2+snf)/fP3;
+ Double_t xv=(fP2+snf)/fP4;
return PropagateTo(xv,x0,rho,pm);
}
Double_t k40=fC40*r00+fC41*r01, k41=fC40*r01+fC41*r11;
Double_t dy=c->GetY() - fP0, dz=c->GetZ() - fP1;
- Double_t cur=fP3 + k30*dy + k31*dz, eta=fP2 + k20*dy + k21*dz;
+ Double_t cur=fP4 + k40*dy + k41*dz, eta=fP2 + k20*dy + k21*dz;
if (TMath::Abs(cur*fX-eta) >= 0.99999) {
Int_t n=GetNumberOfClusters();
if (n>4) cerr<<n<<" AliTPCtrack warning: Filtering failed !\n";
fP0 += k00*dy + k01*dz;
fP1 += k10*dy + k11*dz;
fP2 = eta;
- fP3 = cur;
- fP4 += k40*dy + k41*dz;
+ fP3 += k30*dy + k31*dz;
+ fP4 = cur;
Double_t c01=fC10, c02=fC20, c03=fC30, c04=fC40;
Double_t c12=fC21, c13=fC31, c14=fC41;
fC42-=k20*c04+k21*c14;
fC33-=k30*c03+k31*c13;
- fC43-=k30*c04+k31*c14;
+ fC43-=k40*c03+k41*c13;
fC44-=k40*c04+k41*c14;
Double_t x1=fX, y1=fP0;
Double_t ca=cos(alpha), sa=sin(alpha);
- Double_t r1=fP3*fX - fP2;
+ Double_t r1=fP4*fX - fP2;
fX = x1*ca + y1*sa;
fP0=-x1*sa + y1*ca;
- fP2=fP2*ca + (fP3*y1 + sqrt(1.- r1*r1))*sa;
+ fP2=fP2*ca + (fP4*y1 + sqrt(1.- r1*r1))*sa;
- Double_t r2=fP3*fX - fP2;
+ Double_t r2=fP4*fX - fP2;
if (TMath::Abs(r2) >= 0.99999) {
Int_t n=GetNumberOfClusters();
if (n>4) cerr<<n<<" AliTPCtrack warning: Rotation failed !\n";
return 0;
}
- Double_t y0=fP0 + sqrt(1.- r2*r2)/fP3;
- if ((fP0-y0)*fP3 >= 0.) {
+ Double_t y0=fP0 + sqrt(1.- r2*r2)/fP4;
+ if ((fP0-y0)*fP4 >= 0.) {
Int_t n=GetNumberOfClusters();
if (n>4) cerr<<n<<" AliTPCtrack warning: Rotation failed !!!\n";
return 0;
}
//f = F - 1
- Double_t f00=ca-1, f23=(y1 - r1*x1/sqrt(1.- r1*r1))*sa,
- f20=fP3*sa, f22=(ca + sa*r1/sqrt(1.- r1*r1))-1;
+ Double_t f00=ca-1, f24=(y1 - r1*x1/sqrt(1.- r1*r1))*sa,
+ f20=fP4*sa, f22=(ca + sa*r1/sqrt(1.- r1*r1))-1;
//b = C*ft
- Double_t b00=fC00*f00, b02=fC00*f20+fC30*f23+fC20*f22;
- Double_t b10=fC10*f00, b12=fC10*f20+fC31*f23+fC21*f22;
- Double_t b20=fC20*f00, b22=fC20*f20+fC32*f23+fC22*f22;
- Double_t b30=fC30*f00, b32=fC30*f20+fC33*f23+fC32*f22;
- Double_t b40=fC40*f00, b42=fC40*f20+fC43*f23+fC42*f22;
+ Double_t b00=fC00*f00, b02=fC00*f20+fC40*f24+fC20*f22;
+ Double_t b10=fC10*f00, b12=fC10*f20+fC41*f24+fC21*f22;
+ Double_t b20=fC20*f00, b22=fC20*f20+fC42*f24+fC22*f22;
+ Double_t b30=fC30*f00, b32=fC30*f20+fC43*f24+fC32*f22;
+ Double_t b40=fC40*f00, b42=fC40*f20+fC44*f24+fC42*f22;
//a = f*b = f*C*ft
- Double_t a00=f00*b00, a02=f00*b02, a22=f20*b02+f23*b32+f22*b22;
+ Double_t a00=f00*b00, a02=f00*b02, a22=f20*b02+f24*b42+f22*b22;
// *** Double_t dy2=fCyy;
return 1;
}
-//_____________________________________________________________________________
-void AliTPCtrack::CookLabel(AliTPCClustersArray *ca) {
- //-----------------------------------------------------------------
- // This function cooks the track label. If label<0, this track is fake.
- //-----------------------------------------------------------------
- Int_t n=GetNumberOfClusters();
- Int_t *lb=new Int_t[n];
- Int_t *mx=new Int_t[n];
- AliTPCcluster **clusters=new AliTPCcluster*[n];
-
- Int_t i;
- Int_t sec,row,ncl;
- for (i=0; i<n; i++) {
- lb[i]=mx[i]=0;
- GetCluster(i,sec,row,ncl);
- AliTPCClustersRow *clrow=ca->GetRow(sec,row);
- clusters[i]=(AliTPCcluster*)(*clrow)[ncl];
- }
-
- Int_t lab=123456789;
- for (i=0; i<n; i++) {
- AliTPCcluster *c=clusters[i];
- lab=TMath::Abs(c->GetLabel(0));
- Int_t j;
- for (j=0; j<n; j++)
- if (lb[j]==lab || mx[j]==0) break;
- lb[j]=lab;
- (mx[j])++;
- }
-
- Int_t max=0;
- for (i=0; i<n; i++)
- if (mx[i]>max) {max=mx[i]; lab=lb[i];}
-
- for (i=0; i<n; i++) {
- AliTPCcluster *c=clusters[i];
- if (TMath::Abs(c->GetLabel(1)) == lab ||
- TMath::Abs(c->GetLabel(2)) == lab ) max++;
- }
-
- SetLabel(-lab);
- if (1.-Float_t(max)/n <= 0.10) {
- //Int_t tail=Int_t(0.08*fN);
- Int_t tail=14;
- max=0;
- for (i=1; i<=tail; i++) {
- AliTPCcluster *c=clusters[n-i];
- if (lab == TMath::Abs(c->GetLabel(0)) ||
- lab == TMath::Abs(c->GetLabel(1)) ||
- lab == TMath::Abs(c->GetLabel(2))) max++;
- }
- if (max >= Int_t(0.5*tail)) SetLabel(lab);
- }
-/*
-if (lab==111)
- for (i=0; i<n; i++) {
- AliTPCcluster *c=clusters[i];
- cerr<<i<<' '<<c->GetLabel(0)<<endl;
- }
-*/
- delete[] lb;
- delete[] mx;
- delete[] clusters;
+void AliTPCtrack::ResetCovariance() {
+ //------------------------------------------------------------------
+ //This function makes a track forget its history :)
+ //------------------------------------------------------------------
+
+ fC00*=10.;
+ fC10=0.; fC11*=10.;
+ fC20=0.; fC21=0.; fC22*=10.;
+ fC30=0.; fC31=0.; fC32=0.; fC33*=10.;
+ fC40=0.; fC41=0.; fC42=0.; fC43=0.; fC44*=10.;
+
}
#include <AliKalmanTrack.h>
#include <TMath.h>
+#include "AliTPCreco.h"
+
class AliTPCClustersArray;
class AliTPCcluster;
AliTPCtrack():AliKalmanTrack(){}
AliTPCtrack(UInt_t index, const Double_t xx[5],
const Double_t cc[15], Double_t xr, Double_t alpha);
+ AliTPCtrack(const AliKalmanTrack& t, Double_t alpha);
AliTPCtrack(const AliTPCtrack& t);
Int_t PropagateToVertex(
Double_t x0=36.66,Double_t rho=1.2e-3,Double_t pm=0.139);
Int_t Rotate(Double_t angle);
- void CookLabel(AliTPCClustersArray *carray);
void SetdEdx(Float_t dedx) {fdEdx=dedx;}
Double_t GetX() const {return fX;}
Double_t GetY() const {return fP0;}
Double_t GetZ() const {return fP1;}
- Double_t GetSnp() const {return fX*fP3 - fP2;}
- Double_t Get1Pt() const {return fP3*GetConvConst();}
- Double_t GetTgl() const {return fP4;}
+ Double_t GetSnp() const {return fX*fP4 - fP2;}
+ Double_t Get1Pt() const {return fP4*GetConvConst();}
+ Double_t GetTgl() const {return fP3;}
Double_t GetSigmaY2() const {return fC00;}
Double_t GetSigmaZ2() const {return fC11;}
void GetExternalParameters(Double_t& xr, Double_t x[5]) const ;
void GetExternalCovariance(Double_t cov[15]) const ;
+ Int_t GetClusterIndex(Int_t i) const {return fIndex[i];}
+
//******** To be removed next release !!! **************
Double_t GetEta() const {return fP2;}
- Double_t GetC() const {return fP3;}
+ Double_t GetC() const {return fP4;}
void GetCovariance(Double_t cc[15]) const {
cc[0 ]=fC00;
cc[1 ]=fC10; cc[2 ]=fC11;
cc[3 ]=fC20; cc[4 ]=fC21; cc[5 ]=fC22;
- cc[6 ]=fC30; cc[7 ]=fC31; cc[8 ]=fC32; cc[9 ]=fC33;
- cc[10]=fC40; cc[11]=fC41; cc[12]=fC42; cc[13]=fC43; cc[14]=fC44;
+ cc[6 ]=fC40; cc[7 ]=fC41; cc[8 ]=fC42; cc[9 ]=fC44;
+ cc[10]=fC30; cc[11]=fC31; cc[12]=fC32; cc[13]=fC43; cc[14]=fC33;
}
//******************************************************
- void GetCluster(Int_t i, Int_t &sec, Int_t &row, Int_t &ncl) const;
-
virtual Double_t GetPredictedChi2(const AliCluster *cluster) const;
Int_t PropagateTo(Double_t xr,
Double_t x0=28.94,Double_t rho=0.9e-3,Double_t pm=0.139);
Int_t Update(const AliCluster* c, Double_t chi2, UInt_t i);
+ void ResetCovariance();
private:
Double_t fX; // X-coordinate of this track (reference plane)
Double_t fP0; // Y-coordinate of a track
Double_t fP1; // Z-coordinate of a track
Double_t fP2; // C*x0
- Double_t fP3; // track curvature
- Double_t fP4; // tangent of the track momentum dip angle
+ Double_t fP3; // tangent of the track momentum dip angle
+ Double_t fP4; // track curvature
Double_t fC00; // covariance
Double_t fC10, fC11; // matrix
Double_t fC30, fC31, fC32, fC33; // track
Double_t fC40, fC41, fC42, fC43, fC44; // parameters
- UInt_t fIndex[200]; // indices of associated clusters
+ UInt_t fIndex[kMaxRow]; // indices of associated clusters
ClassDef(AliTPCtrack,1) // Time Projection Chamber reconstructed tracks
};
-
-inline
-void AliTPCtrack::GetCluster(Int_t i,Int_t &sec,Int_t &row,Int_t &ncl) const {
- //return sector, pad row and the index of the i-th cluster of this track
- Int_t index=fIndex[i];
- sec=(index&0xff000000)>>24;
- row=(index&0x00ff0000)>>16;
- ncl=(index&0x0000ffff)>>00;
-}
-
inline
void AliTPCtrack::GetExternalParameters(Double_t& xr, Double_t x[5]) const {
//---------------------------------------------------------------------
/*
$Log$
-Revision 1.6 2001/03/13 14:25:47 hristov
-New design of tracking classes (Yu.Belikov)
-
Revision 1.5 2000/12/20 07:51:59 kowal2
Changes suggested by Alessandra and Paolo to avoid overlapped
data fields in encapsulated classes.
#include "AliTPCtracker.h"
#include "AliTPCcluster.h"
-#include "AliTPCClustersArray.h"
+#include "AliTPCParam.h"
#include "AliTPCClustersRow.h"
-const AliTPCParam *AliTPCtracker::AliTPCSector::fgParam;
+//_____________________________________________________________________________
+AliTPCtracker::AliTPCtracker(const AliTPCParam *par):
+fkNIS(par->GetNInnerSector()/2),
+fkNOS(par->GetNOuterSector()/2)
+{
+ //---------------------------------------------------------------------
+ // The main TPC tracker constructor
+ //---------------------------------------------------------------------
+ fInnerSec=new AliTPCSector[fkNIS];
+ fOuterSec=new AliTPCSector[fkNOS];
+
+ Int_t i;
+ for (i=0; i<fkNIS; i++) fInnerSec[i].Setup(par,0);
+ for (i=0; i<fkNOS; i++) fOuterSec[i].Setup(par,1);
+
+ fN=0; fSectors=0;
+
+ fClustersArray.Setup(par);
+ fClustersArray.SetClusterType("AliTPCcluster");
+ fClustersArray.ConnectTree("Segment Tree");
+
+ fSeeds=0;
+}
+
+//_____________________________________________________________________________
+AliTPCtracker::~AliTPCtracker() {
+ //------------------------------------------------------------------
+ // TPC tracker destructor
+ //------------------------------------------------------------------
+ delete[] fInnerSec;
+ delete[] fOuterSec;
+ delete fSeeds;
+}
//_____________________________________________________________________________
Double_t SigmaY2(Double_t r, Double_t tgl, Double_t pt)
}
//_____________________________________________________________________________
-Int_t AliTPCtracker::FindProlongation(AliTPCseed& t, const AliTPCSector *sec,
-Int_t s, Int_t rf)
-{
+void AliTPCtracker::LoadOuterSectors() {
+ //-----------------------------------------------------------------
+ // This function fills outer TPC sectors with clusters.
+ //-----------------------------------------------------------------
+ UInt_t index;
+ Int_t j=Int_t(fClustersArray.GetTree()->GetEntries());
+ for (Int_t i=0; i<j; i++) {
+ AliSegmentID *s=fClustersArray.LoadEntry(i);
+ Int_t sec,row;
+ AliTPCParam *par=(AliTPCParam*)fClustersArray.GetParam();
+ par->AdjustSectorRow(s->GetID(),sec,row);
+ if (sec<fkNIS*2) continue;
+ AliTPCClustersRow *clrow=fClustersArray.GetRow(sec,row);
+ Int_t ncl=clrow->GetArray()->GetEntriesFast();
+ while (ncl--) {
+ AliTPCcluster *c=(AliTPCcluster*)(*clrow)[ncl];
+ index=(((sec<<8)+row)<<16)+ncl;
+ fOuterSec[(sec-fkNIS*2)%fkNOS][row].InsertCluster(c,index);
+ }
+ }
+
+ fN=fkNOS;
+ fSectors=fOuterSec;
+}
+
+//_____________________________________________________________________________
+void AliTPCtracker::UnloadOuterSectors() {
+ //-----------------------------------------------------------------
+ // This function clears outer TPC sectors.
+ //-----------------------------------------------------------------
+ Int_t nup=fOuterSec->GetNRows();
+ for (Int_t i=0; i<fkNOS; i++) {
+ for (Int_t j=0; j<nup; j++) {
+ if (fClustersArray.GetRow(i+fkNIS*2,j))
+ fClustersArray.ClearRow(i+fkNIS*2,j);
+ if (fClustersArray.GetRow(i+fkNIS*2+fkNOS,j))
+ fClustersArray.ClearRow(i+fkNIS*2+fkNOS,j);
+ }
+ }
+
+ fN=0;
+ fSectors=0;
+}
+
+//_____________________________________________________________________________
+void AliTPCtracker::LoadInnerSectors() {
+ //-----------------------------------------------------------------
+ // This function fills inner TPC sectors with clusters.
+ //-----------------------------------------------------------------
+ UInt_t index;
+ Int_t j=Int_t(fClustersArray.GetTree()->GetEntries());
+ for (Int_t i=0; i<j; i++) {
+ AliSegmentID *s=fClustersArray.LoadEntry(i);
+ Int_t sec,row;
+ AliTPCParam *par=(AliTPCParam*)fClustersArray.GetParam();
+ par->AdjustSectorRow(s->GetID(),sec,row);
+ if (sec>=fkNIS*2) continue;
+ AliTPCClustersRow *clrow=fClustersArray.GetRow(sec,row);
+ Int_t ncl=clrow->GetArray()->GetEntriesFast();
+ while (ncl--) {
+ AliTPCcluster *c=(AliTPCcluster*)(*clrow)[ncl];
+ index=(((sec<<8)+row)<<16)+ncl;
+ fInnerSec[sec%fkNIS][row].InsertCluster(c,index);
+ }
+ }
+
+ fN=fkNIS;
+ fSectors=fInnerSec;
+}
+
+//_____________________________________________________________________________
+void AliTPCtracker::UnloadInnerSectors() {
+ //-----------------------------------------------------------------
+ // This function clears inner TPC sectors.
+ //-----------------------------------------------------------------
+ Int_t nlow=fInnerSec->GetNRows();
+ for (Int_t i=0; i<fkNIS; i++) {
+ for (Int_t j=0; j<nlow; j++) {
+ if (fClustersArray.GetRow(i,j)) fClustersArray.ClearRow(i,j);
+ if (fClustersArray.GetRow(i+fkNIS,j)) fClustersArray.ClearRow(i+fkNIS,j);
+ }
+ }
+
+ fN=0;
+ fSectors=0;
+}
+
+//_____________________________________________________________________________
+Int_t AliTPCtracker::FollowProlongation(AliTPCseed& t, Int_t rf) {
//-----------------------------------------------------------------
// This function tries to find a track prolongation.
//-----------------------------------------------------------------
- const Int_t kSKIP=(t.GetNumberOfClusters()<10) ? 10 :
- Int_t(0.5*sec->GetNRows());
+ Double_t xt=t.GetX();
+ const Int_t kSKIP=(t.GetNumberOfClusters()<10) ? kRowsToSkip :
+ Int_t(0.5*fSectors->GetNRows());
Int_t tryAgain=kSKIP;
- Double_t alpha=sec->GetAlpha();
- Int_t ns=Int_t(2*TMath::Pi()/alpha+0.5);
- for (Int_t nr=sec->GetRowNumber(t.GetX())-1; nr>=rf; nr--) {
- Double_t x=sec->GetX(nr), ymax=sec->GetMaxY(nr);
+ Double_t alpha=t.GetAlpha() - fSectors->GetAlphaShift();
+ if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
+ if (alpha < 0. ) alpha += 2.*TMath::Pi();
+ Int_t s=Int_t(alpha/fSectors->GetAlpha())%fN;
+
+ for (Int_t nr=fSectors->GetRowNumber(xt)-1; nr>=rf; nr--) {
+ Double_t x=fSectors->GetX(nr), ymax=fSectors->GetMaxY(nr);
if (!t.PropagateTo(x)) return 0;
AliTPCcluster *cl=0;
UInt_t index=0;
- Double_t maxchi2=12.;
- const AliTPCRow &krow=sec[s][nr];
- Double_t sy2=SigmaY2(t.GetX(),t.GetTgl(),1./t.Get1Pt());
+ Double_t maxchi2=kMaxCHI2;
+ const AliTPCRow &krow=fSectors[s][nr];
+ Double_t pt=t.GetConvConst()/(100/0.299792458/0.2)/t.Get1Pt();
+ Double_t sy2=SigmaY2(t.GetX(),t.GetTgl(),pt);
Double_t sz2=SigmaZ2(t.GetX(),t.GetTgl());
Double_t road=4.*sqrt(t.GetSigmaY2() + sy2), y=t.GetY(), z=t.GetZ();
- if (road>30) {
+ if (road>kMaxROAD) {
if (t.GetNumberOfClusters()>4)
cerr<<t.GetNumberOfClusters()
<<"FindProlongation warning: Too broad road !\n";
}
}
if (cl) {
- Float_t l=sec->GetPadPitchWidth();
+ Float_t l=fSectors->GetPadPitchWidth();
t.SetSampledEdx(cl->GetQ()/l,t.GetNumberOfClusters());
if (!t.Update(cl,maxchi2,index)) {
if (!tryAgain--) return 0;
} else {
if (tryAgain==0) break;
if (y > ymax) {
- s = (s+1) % ns;
- if (!t.Rotate(alpha)) return 0;
+ s = (s+1) % fN;
+ if (!t.Rotate(fSectors->GetAlpha())) return 0;
} else if (y <-ymax) {
- s = (s-1+ns) % ns;
- if (!t.Rotate(-alpha)) return 0;
+ s = (s-1+fN) % fN;
+ if (!t.Rotate(-fSectors->GetAlpha())) return 0;
}
tryAgain--;
}
}
return 1;
+}
+
+//_____________________________________________________________________________
+Int_t AliTPCtracker::FollowBackProlongation
+(AliTPCseed& seed, const AliTPCtrack &track) {
+ //-----------------------------------------------------------------
+ // This function propagates tracks back through the TPC
+ //-----------------------------------------------------------------
+ Double_t alpha=seed.GetAlpha() - fSectors->GetAlphaShift();
+ if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
+ if (alpha < 0. ) alpha += 2.*TMath::Pi();
+ Int_t s=Int_t(alpha/fSectors->GetAlpha())%fN;
+
+ Int_t idx=-1, sec=-1, row=-1;
+ Int_t nc=seed.GetLabel(); //index of the cluster to start with
+ if (nc--) {
+ idx=track.GetClusterIndex(nc);
+ sec=(idx&0xff000000)>>24; row=(idx&0x00ff0000)>>16;
+ }
+ if (fSectors==fInnerSec) { if (sec >= 2*fkNIS) row=-1; }
+ else { if (sec < 2*fkNIS) row=-1; }
+
+ Int_t nr=fSectors->GetNRows();
+ for (Int_t i=0; i<nr; i++) {
+ Double_t x=fSectors->GetX(i), ymax=fSectors->GetMaxY(i);
+
+ if (!seed.PropagateTo(x)) return 0;
+
+ Double_t y=seed.GetY();
+ if (y > ymax) {
+ s = (s+1) % fN;
+ if (!seed.Rotate(fSectors->GetAlpha())) return 0;
+ } else if (y <-ymax) {
+ s = (s-1+fN) % fN;
+ if (!seed.Rotate(-fSectors->GetAlpha())) return 0;
+ }
+
+ AliTPCcluster *cl=0;
+ Int_t index=0;
+ Double_t maxchi2=kMaxCHI2;
+ Double_t pt=seed.GetConvConst()/(100/0.299792458/0.2)/seed.Get1Pt();
+ Double_t sy2=SigmaY2(seed.GetX(),seed.GetTgl(),pt);
+ Double_t sz2=SigmaZ2(seed.GetX(),seed.GetTgl());
+ Double_t road=4.*sqrt(seed.GetSigmaY2() + sy2), z=seed.GetZ();
+ if (road>kMaxROAD) {
+ cerr<<seed.GetNumberOfClusters()
+ <<"AliTPCtracker::FollowBackProlongation: Too broad road !\n";
+ return 0;
+ }
+
+
+ Int_t accepted=seed.GetNumberOfClusters();
+ if (row==i) {
+ //try to accept already found cluster
+ AliTPCcluster *c=(AliTPCcluster*)GetCluster(idx);
+ Double_t chi2;
+ if ((chi2=seed.GetPredictedChi2(c))<maxchi2 || accepted<27) {
+ index=idx; cl=c; maxchi2=chi2;
+ } //else cerr<<"AliTPCtracker::FollowBackProlongation: oulier !\n";
+
+ if (nc--) {
+ idx=track.GetClusterIndex(nc);
+ sec=(idx&0xff000000)>>24; row=(idx&0x00ff0000)>>16;
+ }
+ if (fSectors==fInnerSec) { if (sec >= 2*fkNIS) row=-1; }
+ else { if (sec < 2*fkNIS) row=-1; }
+ }
+ if (!cl) {
+ //try to fill the gap
+ const AliTPCRow &krow=fSectors[s][i];
+ if (accepted>27)
+ if (krow) {
+ for (Int_t i=krow.Find(y-road); i<krow; i++) {
+ AliTPCcluster *c=(AliTPCcluster*)(krow[i]);
+ if (c->GetY() > y+road) break;
+ if (c->IsUsed()) continue;
+ if ((c->GetZ()-z)*(c->GetZ()-z)>16.*(seed.GetSigmaZ2()+sz2)) continue;
+ Double_t chi2=seed.GetPredictedChi2(c);
+ if (chi2 > maxchi2) continue;
+ maxchi2=chi2;
+ cl=c;
+ index=krow.GetIndex(i);
+ }
+ }
+ }
+
+ if (cl) {
+ Float_t l=fSectors->GetPadPitchWidth();
+ seed.SetSampledEdx(cl->GetQ()/l,seed.GetNumberOfClusters());
+ seed.Update(cl,maxchi2,index);
+ }
+
+ }
+
+ seed.SetLabel(nc);
+
+ return 1;
}
//_____________________________________________________________________________
-void
-AliTPCtracker::MakeSeeds(TObjArray& seeds,const AliTPCSector *sec, Int_t max,
-Int_t i1, Int_t i2)
-{
+void AliTPCtracker::MakeSeeds(Int_t i1, Int_t i2) {
//-----------------------------------------------------------------
// This function creates track seeds.
//-----------------------------------------------------------------
+ if (fSeeds==0) fSeeds=new TObjArray(15000);
+
Double_t x[5], c[15];
- Double_t alpha=sec->GetAlpha(), shift=sec->GetAlphaShift();
+ Double_t alpha=fOuterSec->GetAlpha(), shift=fOuterSec->GetAlphaShift();
Double_t cs=cos(alpha), sn=sin(alpha);
- Double_t x1 =sec->GetX(i1);
- Double_t xx2=sec->GetX(i2);
+ Double_t x1 =fOuterSec->GetX(i1);
+ Double_t xx2=fOuterSec->GetX(i2);
- for (Int_t ns=0; ns<max; ns++) {
- Int_t nl=sec[(ns-1+max)%max][i2];
- Int_t nm=sec[ns][i2];
- Int_t nu=sec[(ns+1)%max][i2];
- const AliTPCRow& kr1=sec[ns][i1];
+ for (Int_t ns=0; ns<fkNOS; ns++) {
+ Int_t nl=fOuterSec[(ns-1+fkNOS)%fkNOS][i2];
+ Int_t nm=fOuterSec[ns][i2];
+ Int_t nu=fOuterSec[(ns+1)%fkNOS][i2];
+ const AliTPCRow& kr1=fOuterSec[ns][i1];
for (Int_t is=0; is < kr1; is++) {
Double_t y1=kr1[is]->GetY(), z1=kr1[is]->GetZ();
for (Int_t js=0; js < nl+nm+nu; js++) {
Double_t x3=0.,y3=0.;
if (js<nl) {
- const AliTPCRow& kr2=sec[(ns-1+max)%max][i2];
+ const AliTPCRow& kr2=fOuterSec[(ns-1+fkNOS)%fkNOS][i2];
kcl=kr2[js];
y2=kcl->GetY(); z2=kcl->GetZ();
x2= xx2*cs+y2*sn;
y2=-xx2*sn+y2*cs;
} else
if (js<nl+nm) {
- const AliTPCRow& kr2=sec[ns][i2];
+ const AliTPCRow& kr2=fOuterSec[ns][i2];
kcl=kr2[js-nl];
x2=xx2; y2=kcl->GetY(); z2=kcl->GetZ();
} else {
- const AliTPCRow& kr2=sec[(ns+1)%max][i2];
+ const AliTPCRow& kr2=fOuterSec[(ns+1)%fkNOS][i2];
kcl=kr2[js-nl-nm];
y2=kcl->GetY(); z2=kcl->GetZ();
x2=xx2*cs-y2*sn;
x[0]=y1;
x[1]=z1;
- x[3]=f1(x1,y1,x2,y2,x3,y3);
- if (TMath::Abs(x[3]) >= 0.0066) continue;
+ x[4]=f1(x1,y1,x2,y2,x3,y3);
+ if (TMath::Abs(x[4]) >= 0.0066) continue;
x[2]=f2(x1,y1,x2,y2,x3,y3);
- //if (TMath::Abs(x[3]*x1-x[2]) >= 0.99999) continue;
- x[4]=f3(x1,y1,x2,y2,z1,z2);
- if (TMath::Abs(x[4]) > 1.2) continue;
+ //if (TMath::Abs(x[4]*x1-x[2]) >= 0.99999) continue;
+ x[3]=f3(x1,y1,x2,y2,z1,z2);
+ if (TMath::Abs(x[3]) > 1.2) continue;
Double_t a=asin(x[2]);
- Double_t zv=z1 - x[4]/x[3]*(a+asin(x[3]*x1-x[2]));
+ Double_t zv=z1 - x[3]/x[4]*(a+asin(x[4]*x1-x[2]));
if (TMath::Abs(zv)>10.) continue;
Double_t sy1=kr1[is]->GetSigmaY2(), sz1=kr1[is]->GetSigmaZ2();
Double_t sy2=kcl->GetSigmaY2(), sz2=kcl->GetSigmaZ2();
Double_t sy3=100*0.025, sy=0.1, sz=0.1;
- Double_t f30=(f1(x1,y1+sy,x2,y2,x3,y3)-x[3])/sy;
- Double_t f32=(f1(x1,y1,x2,y2+sy,x3,y3)-x[3])/sy;
- Double_t f34=(f1(x1,y1,x2,y2,x3,y3+sy)-x[3])/sy;
+ Double_t f40=(f1(x1,y1+sy,x2,y2,x3,y3)-x[4])/sy;
+ Double_t f42=(f1(x1,y1,x2,y2+sy,x3,y3)-x[4])/sy;
+ Double_t f43=(f1(x1,y1,x2,y2,x3,y3+sy)-x[4])/sy;
Double_t f20=(f2(x1,y1+sy,x2,y2,x3,y3)-x[2])/sy;
Double_t f22=(f2(x1,y1,x2,y2+sy,x3,y3)-x[2])/sy;
- Double_t f24=(f2(x1,y1,x2,y2,x3,y3+sy)-x[2])/sy;
- Double_t f40=(f3(x1,y1+sy,x2,y2,z1,z2)-x[4])/sy;
- Double_t f41=(f3(x1,y1,x2,y2,z1+sz,z2)-x[4])/sz;
- Double_t f42=(f3(x1,y1,x2,y2+sy,z1,z2)-x[4])/sy;
- Double_t f43=(f3(x1,y1,x2,y2,z1,z2+sz)-x[4])/sz;
+ Double_t f23=(f2(x1,y1,x2,y2,x3,y3+sy)-x[2])/sy;
+ Double_t f30=(f3(x1,y1+sy,x2,y2,z1,z2)-x[3])/sy;
+ Double_t f31=(f3(x1,y1,x2,y2,z1+sz,z2)-x[3])/sz;
+ Double_t f32=(f3(x1,y1,x2,y2+sy,z1,z2)-x[3])/sy;
+ Double_t f34=(f3(x1,y1,x2,y2,z1,z2+sz)-x[3])/sz;
c[0]=sy1;
c[1]=0.; c[2]=sz1;
- c[3]=f20*sy1; c[4]=0.; c[5]=f20*sy1*f20+f22*sy2*f22+f24*sy3*f24;
- c[6]=f30*sy1; c[7]=0.; c[8]=f30*sy1*f20+f32*sy2*f22+f34*sy3*f24;
- c[9]=f30*sy1*f30+f32*sy2*f32+f34*sy3*f34;
- c[10]=f40*sy1; c[11]=f41*sz1; c[12]=f40*sy1*f20+f42*sy2*f22;
- c[13]=f40*sy1*f30+f42*sy2*f32;
- c[14]=f40*sy1*f40+f41*sz1*f41+f42*sy2*f42+f43*sz2*f43;
+ c[3]=f20*sy1; c[4]=0.; c[5]=f20*sy1*f20+f22*sy2*f22+f23*sy3*f23;
+ c[6]=f30*sy1; c[7]=f31*sz1; c[8]=f30*sy1*f20+f32*sy2*f22;
+ c[9]=f30*sy1*f30+f31*sz1*f31+f32*sy2*f32+f34*sz2*f34;
+ c[10]=f40*sy1; c[11]=0.; c[12]=f40*sy1*f20+f42*sy2*f22+f43*sy3*f23;
+ c[13]=f30*sy1*f40+f32*sy2*f42;
+ c[14]=f40*sy1*f40+f42*sy2*f42+f43*sy3*f43;
UInt_t index=kr1.GetIndex(is);
AliTPCseed *track=new AliTPCseed(index, x, c, x1, ns*alpha+shift);
- Float_t l=sec->GetPadPitchWidth();
+ Float_t l=fOuterSec->GetPadPitchWidth();
track->SetSampledEdx(kr1[is]->GetQ()/l,0);
- Int_t rc=FindProlongation(*track,sec,ns,i2);
+ Int_t rc=FollowProlongation(*track, i2);
if (rc==0 || track->GetNumberOfClusters()<(i1-i2)/2) delete track;
- else seeds.AddLast(track);
+ else fSeeds->AddLast(track);
}
}
}
}
//_____________________________________________________________________________
-Int_t AliTPCtracker::Clusters2Tracks(const AliTPCParam *par, TFile *of) {
+Int_t AliTPCtracker::ReadSeeds(const TFile *inp) {
//-----------------------------------------------------------------
- // This is a track finder.
+ // This function reades track seeds.
//-----------------------------------------------------------------
TDirectory *savedir=gDirectory;
- if (!of->IsOpen()) {
- cerr<<"AliTPCtracker::Clusters2Tracks(): output file not open !\n";
+ TFile *in=(TFile*)inp;
+ if (!in->IsOpen()) {
+ cerr<<"AliTPCtracker::ReadSeeds(): input file is not open !\n";
return 1;
}
- AliTPCClustersArray carray;
- carray.Setup(par);
- carray.SetClusterType("AliTPCcluster");
- carray.ConnectTree("Segment Tree");
+ in->cd();
+ TTree *seedTree=(TTree*)in->Get("Seeds");
+ if (!seedTree) {
+ cerr<<"AliTPCtracker::ReadSeeds(): ";
+ cerr<<"can't get a tree with track seeds !\n";
+ return 2;
+ }
+ AliTPCtrack *seed=new AliTPCtrack;
+ seedTree->SetBranchAddress("tracks",&seed);
+
+ if (fSeeds==0) fSeeds=new TObjArray(15000);
- of->cd();
- TTree tracktree("TreeT","Tree with TPC tracks");
- AliTPCtrack *iotrack=0;
- tracktree.Branch("tracks","AliTPCtrack",&iotrack,32000,0);
+ Int_t n=(Int_t)seedTree->GetEntries();
+ for (Int_t i=0; i<n; i++) {
+ seedTree->GetEvent(i);
+ fSeeds->AddLast(new AliTPCseed(*seed,seed->GetAlpha()));
+ }
+
+ delete seed;
+ savedir->cd();
- AliTPCSector::SetParam(par);
+ return 0;
+}
- const Int_t kNIS=par->GetNInnerSector()/2;
- AliTPCSSector *ssec=new AliTPCSSector[kNIS];
- Int_t nlow=ssec->GetNRows();
+//_____________________________________________________________________________
+Int_t AliTPCtracker::Clusters2Tracks(const TFile *inp, TFile *out) {
+ //-----------------------------------------------------------------
+ // This is a track finder.
+ //-----------------------------------------------------------------
+ TDirectory *savedir=gDirectory;
- const Int_t kNOS=par->GetNOuterSector()/2;
- AliTPCLSector *lsec=new AliTPCLSector[kNOS];
- Int_t nup=lsec->GetNRows();
-
- //Load outer sectors
- UInt_t index;
- Int_t i,j;
+ if (inp) {
+ TFile *in=(TFile*)inp;
+ if (!in->IsOpen()) {
+ cerr<<"AliTPCtracker::Clusters2Tracks(): input file is not open !\n";
+ return 1;
+ }
+ }
- j=Int_t(carray.GetTree()->GetEntries());
- for (i=0; i<j; i++) {
- AliSegmentID *s=carray.LoadEntry(i);
- Int_t sec,row;
- par->AdjustSectorRow(s->GetID(),sec,row);
- if (sec<kNIS*2) continue;
- AliTPCClustersRow *clrow=carray.GetRow(sec,row);
- Int_t ncl=clrow->GetArray()->GetEntriesFast();
- while (ncl--) {
- AliTPCcluster *c=(AliTPCcluster*)(*clrow)[ncl];
- index=(((sec<<8)+row)<<16)+ncl;
- lsec[(sec-kNIS*2)%kNOS][row].InsertCluster(c,index);
- }
+ if (!out->IsOpen()) {
+ cerr<<"AliTPCtracker::Clusters2Tracks(): output file is not open !\n";
+ return 2;
}
+ out->cd();
+ TTree tracktree("TPCf","Tree with TPC tracks");
+ AliTPCtrack *iotrack=0;
+ tracktree.Branch("tracks","AliTPCtrack",&iotrack,32000,0);
+
+ LoadOuterSectors();
+
//find track seeds
- TObjArray seeds(20000);
+ Int_t nup=fOuterSec->GetNRows(), nlow=fInnerSec->GetNRows();
Int_t nrows=nlow+nup;
- Int_t gap=Int_t(0.125*nrows), shift=Int_t(0.5*gap);
- MakeSeeds(seeds, lsec, kNOS, nup-1, nup-1-gap);
- MakeSeeds(seeds, lsec, kNOS, nup-1-shift, nup-1-shift-gap);
- seeds.Sort();
+ if (fSeeds==0) {
+ Int_t gap=Int_t(0.125*nrows), shift=Int_t(0.5*gap);
+ MakeSeeds(nup-1, nup-1-gap);
+ MakeSeeds(nup-1-shift, nup-1-shift-gap);
+ }
+ fSeeds->Sort();
//tracking in outer sectors
- Int_t nseed=seeds.GetEntriesFast();
+ Int_t nseed=fSeeds->GetEntriesFast();
+ Int_t i;
for (i=0; i<nseed; i++) {
- AliTPCseed *pt=(AliTPCseed*)seeds.UncheckedAt(i), &t=*pt;
- Double_t alpha=t.GetAlpha();
- if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
- if (alpha < 0. ) alpha += 2.*TMath::Pi();
- Int_t ns=Int_t(alpha/lsec->GetAlpha())%kNOS;
-
- if (FindProlongation(t,lsec,ns)) {
- t.UseClusters(&carray);
+ AliTPCseed *pt=(AliTPCseed*)fSeeds->UncheckedAt(i), &t=*pt;
+ if (FollowProlongation(t)) {
+ UseClusters(&t);
continue;
}
- delete seeds.RemoveAt(i);
+ delete fSeeds->RemoveAt(i);
}
-
- //unload outer sectors
- for (i=0; i<kNOS; i++) {
- for (j=0; j<nup; j++) {
- if (carray.GetRow(i+kNIS*2,j)) carray.ClearRow(i+kNIS*2,j);
- if (carray.GetRow(i+kNIS*2+kNOS,j)) carray.ClearRow(i+kNIS*2+kNOS,j);
- }
- }
-
- //load inner sectors
- j=Int_t(carray.GetTree()->GetEntries());
- for (i=0; i<j; i++) {
- AliSegmentID *s=carray.LoadEntry(i);
- Int_t sec,row;
- par->AdjustSectorRow(s->GetID(),sec,row);
- if (sec>=kNIS*2) continue;
- AliTPCClustersRow *clrow=carray.GetRow(sec,row);
- Int_t ncl=clrow->GetArray()->GetEntriesFast();
- while (ncl--) {
- AliTPCcluster *c=(AliTPCcluster*)(*clrow)[ncl];
- index=(((sec<<8)+row)<<16)+ncl;
- ssec[sec%kNIS][row].InsertCluster(c,index);
- }
- }
+ UnloadOuterSectors();
//tracking in inner sectors
+ LoadInnerSectors();
Int_t found=0;
for (i=0; i<nseed; i++) {
- AliTPCseed *pt=(AliTPCseed*)seeds.UncheckedAt(i), &t=*pt;
+ AliTPCseed *pt=(AliTPCseed*)fSeeds->UncheckedAt(i), &t=*pt;
if (!pt) continue;
Int_t nc=t.GetNumberOfClusters();
- Double_t alpha=t.GetAlpha() - ssec->GetAlphaShift();
+ Double_t alpha=t.GetAlpha() - fInnerSec->GetAlphaShift();
if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
if (alpha < 0. ) alpha += 2.*TMath::Pi();
- Int_t ns=Int_t(alpha/ssec->GetAlpha())%kNIS;
+ Int_t ns=Int_t(alpha/fInnerSec->GetAlpha())%fkNIS;
- alpha=ns*ssec->GetAlpha() + ssec->GetAlphaShift() - t.GetAlpha();
+ alpha=ns*fInnerSec->GetAlpha() + fInnerSec->GetAlphaShift() - t.GetAlpha();
if (t.Rotate(alpha)) {
- if (FindProlongation(t,ssec,ns)) {
+ if (FollowProlongation(t)) {
if (t.GetNumberOfClusters() >= Int_t(0.4*nrows)) {
t.CookdEdx();
- //t.CookLabel(&carray);
iotrack=pt;
tracktree.Fill();
- t.UseClusters(&carray,nc);
- found++;
+ UseClusters(&t,nc);
+ cerr<<found++<<'\r';
}
}
}
- delete seeds.RemoveAt(i);
+ delete fSeeds->RemoveAt(i);
}
+ UnloadInnerSectors();
+
tracktree.Write();
- //unload inner sectors
- for (i=0; i<kNIS; i++) {
- for (j=0; j<nlow; j++) {
- if (carray.GetRow(i,j)) carray.ClearRow(i,j);
- if (carray.GetRow(i+kNIS,j)) carray.ClearRow(i+kNIS,j);
- }
+ cerr<<"Number of found tracks : "<<found<<endl;
+
+ savedir->cd();
+
+ return 0;
+}
+
+//_____________________________________________________________________________
+Int_t AliTPCtracker::PropagateBack(const TFile *inp, TFile *out) {
+ //-----------------------------------------------------------------
+ // This function propagates tracks back through the TPC.
+ //-----------------------------------------------------------------
+ fSeeds=new TObjArray(15000);
+
+ TFile *in=(TFile*)inp;
+ TDirectory *savedir=gDirectory;
+
+ if (!in->IsOpen()) {
+ cerr<<"AliTPCtracker::PropagateBack(): ";
+ cerr<<"file with back propagated ITS tracks is not open !\n";
+ return 1;
}
- cerr<<"Number of found tracks : "<<found<<endl;
+ if (!out->IsOpen()) {
+ cerr<<"AliTPCtracker::PropagateBack(): ";
+ cerr<<"file for back propagated TPC tracks is not open !\n";
+ return 2;
+ }
+
+ in->cd();
+ TTree *bckTree=(TTree*)in->Get("ITSb");
+ if (!bckTree) {
+ cerr<<"AliTPCtracker::PropagateBack() ";
+ cerr<<"can't get a tree with back propagated ITS tracks !\n";
+ return 3;
+ }
+ AliTPCtrack *bckTrack=new AliTPCtrack;
+ bckTree->SetBranchAddress("tracks",&bckTrack);
+
+ TTree *tpcTree=(TTree*)in->Get("TPCf");
+ if (!tpcTree) {
+ cerr<<"AliTPCtracker::PropagateBack() ";
+ cerr<<"can't get a tree with TPC tracks !\n";
+ return 4;
+ }
+ AliTPCtrack *tpcTrack=new AliTPCtrack;
+ tpcTree->SetBranchAddress("tracks",&tpcTrack);
+
+//*** Prepare an array of tracks to be back propagated
+ Int_t nup=fOuterSec->GetNRows(), nlow=fInnerSec->GetNRows();
+ Int_t nrows=nlow+nup;
+
+ TObjArray tracks(15000);
+ Int_t i=0,j=0;
+ Int_t tpcN=(Int_t)tpcTree->GetEntries();
+ Int_t bckN=(Int_t)bckTree->GetEntries();
+ if (j<bckN) bckTree->GetEvent(j++);
+ for (i=0; i<tpcN; i++) {
+ tpcTree->GetEvent(i);
+ Double_t alpha=tpcTrack->GetAlpha();
+ if (j<bckN &&
+ TMath::Abs(tpcTrack->GetLabel())==TMath::Abs(bckTrack->GetLabel())) {
+ if (!bckTrack->Rotate(alpha-bckTrack->GetAlpha())) continue;
+ fSeeds->AddLast(new AliTPCseed(*bckTrack,bckTrack->GetAlpha()));
+ bckTree->GetEvent(j++);
+ } else {
+ tpcTrack->ResetCovariance();
+ fSeeds->AddLast(new AliTPCseed(*tpcTrack,alpha));
+ }
+ tracks.AddLast(new AliTPCtrack(*tpcTrack));
+ }
+
+ out->cd();
+ TTree backTree("TPCb","Tree with back propagated TPC tracks");
+ AliTPCtrack *otrack=0;
+ backTree.Branch("tracks","AliTPCtrack",&otrack,32000,0);
+
+//*** Back propagation through inner sectors
+ LoadInnerSectors();
+ Int_t nseed=fSeeds->GetEntriesFast();
+ for (i=0; i<nseed; i++) {
+ AliTPCseed *ps=(AliTPCseed*)fSeeds->UncheckedAt(i), &s=*ps;
+ const AliTPCtrack *pt=(AliTPCtrack*)tracks.UncheckedAt(i), &t=*pt;
- delete[] ssec;
- delete[] lsec;
+ Int_t nc=t.GetNumberOfClusters();
+ s.SetLabel(nc-1); //set number of the cluster to start with
+ if (FollowBackProlongation(s,t)) {
+ UseClusters(&s);
+ continue;
+ }
+ delete fSeeds->RemoveAt(i);
+ }
+ UnloadInnerSectors();
+
+//*** Back propagation through outer sectors
+ LoadOuterSectors();
+ Int_t found=0;
+ for (i=0; i<nseed; i++) {
+ AliTPCseed *ps=(AliTPCseed*)fSeeds->UncheckedAt(i), &s=*ps;
+ if (!ps) continue;
+ Int_t nc=s.GetNumberOfClusters();
+ const AliTPCtrack *pt=(AliTPCtrack*)tracks.UncheckedAt(i), &t=*pt;
+
+ Double_t alpha=s.GetAlpha() - fSectors->GetAlphaShift();
+ if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi();
+ if (alpha < 0. ) alpha += 2.*TMath::Pi();
+ Int_t ns=Int_t(alpha/fSectors->GetAlpha())%fN;
+
+ alpha =ns*fSectors->GetAlpha() + fSectors->GetAlphaShift();
+ alpha-=s.GetAlpha();
+
+ if (s.Rotate(alpha)) {
+ if (FollowBackProlongation(s,t)) {
+ if (s.GetNumberOfClusters() >= Int_t(0.4*nrows)) {
+ s.CookdEdx();
+ s.SetLabel(t.GetLabel());
+ UseClusters(&s,nc);
+ otrack=ps;
+ backTree.Fill();
+ cerr<<found++<<'\r';
+ continue;
+ }
+ }
+ }
+ delete fSeeds->RemoveAt(i);
+ }
+ UnloadOuterSectors();
+
+ backTree.Write();
savedir->cd();
+ cerr<<"Number of seeds: "<<nseed<<endl;
+ cerr<<"Number of back propagated ITS tracks: "<<bckN<<endl;
+ cerr<<"Number of back propagated TPC tracks: "<<found<<endl;
+
+ delete bckTrack;
+ delete tpcTrack;
return 0;
}
+//_________________________________________________________________________
+AliCluster *AliTPCtracker::GetCluster(Int_t index) const {
+ //--------------------------------------------------------------------
+ // Return pointer to a given cluster
+ //--------------------------------------------------------------------
+ Int_t sec=(index&0xff000000)>>24;
+ Int_t row=(index&0x00ff0000)>>16;
+ Int_t ncl=(index&0x0000ffff)>>00;
+
+ AliTPCClustersRow *clrow=fClustersArray.GetRow(sec,row);
+ return (AliCluster*)(*clrow)[ncl];
+}
+
+//__________________________________________________________________________
+void AliTPCtracker::CookLabel(AliKalmanTrack *t, Float_t wrong) const {
+ //--------------------------------------------------------------------
+ //This function "cooks" a track label. If label<0, this track is fake.
+ //--------------------------------------------------------------------
+ Int_t noc=t->GetNumberOfClusters();
+ Int_t *lb=new Int_t[noc];
+ Int_t *mx=new Int_t[noc];
+ AliCluster **clusters=new AliCluster*[noc];
+
+ Int_t i;
+ for (i=0; i<noc; i++) {
+ lb[i]=mx[i]=0;
+ Int_t index=t->GetClusterIndex(i);
+ clusters[i]=GetCluster(index);
+ }
+
+ Int_t lab=123456789;
+ for (i=0; i<noc; i++) {
+ AliCluster *c=clusters[i];
+ lab=TMath::Abs(c->GetLabel(0));
+ Int_t j;
+ for (j=0; j<noc; j++) if (lb[j]==lab || mx[j]==0) break;
+ lb[j]=lab;
+ (mx[j])++;
+ }
+
+ Int_t max=0;
+ for (i=0; i<noc; i++) if (mx[i]>max) {max=mx[i]; lab=lb[i];}
+
+ for (i=0; i<noc; i++) {
+ AliCluster *c=clusters[i];
+ if (TMath::Abs(c->GetLabel(1)) == lab ||
+ TMath::Abs(c->GetLabel(2)) == lab ) max++;
+ }
+
+ if ((1.- Float_t(max)/noc) > wrong) lab=-lab;
+
+ else {
+ Int_t tail=Int_t(0.10*noc);
+ max=0;
+ for (i=1; i<=tail; i++) {
+ AliCluster *c=clusters[noc-i];
+ if (lab == TMath::Abs(c->GetLabel(0)) ||
+ lab == TMath::Abs(c->GetLabel(1)) ||
+ lab == TMath::Abs(c->GetLabel(2))) max++;
+ }
+ if (max < Int_t(0.5*tail)) lab=-lab;
+ }
+
+ t->SetLabel(lab);
+
+ delete[] lb;
+ delete[] mx;
+ delete[] clusters;
+}
+
+//_________________________________________________________________________
+void AliTPCtracker::AliTPCSector::Setup(const AliTPCParam *par, Int_t f) {
+ //-----------------------------------------------------------------------
+ // Setup inner sector
+ //-----------------------------------------------------------------------
+ if (f==0) {
+ fAlpha=par->GetInnerAngle();
+ fAlphaShift=par->GetInnerAngleShift();
+ fPadPitchWidth=par->GetInnerPadPitchWidth();
+ fPadPitchLength=par->GetInnerPadPitchLength();
+
+ fN=par->GetNRowLow();
+ fRow=new AliTPCRow[fN];
+ for (Int_t i=0; i<fN; i++) fRow[i].SetX(par->GetPadRowRadiiLow(i));
+ } else {
+ fAlpha=par->GetOuterAngle();
+ fAlphaShift=par->GetOuterAngleShift();
+ fPadPitchWidth=par->GetOuterPadPitchWidth();
+ fPadPitchLength=par->GetOuterPadPitchLength();
+
+ fN=par->GetNRowUp();
+ fRow=new AliTPCRow[fN];
+ for (Int_t i=0; i<fN; i++) fRow[i].SetX(par->GetPadRowRadiiUp(i));
+ }
+}
+
//_________________________________________________________________________
void
AliTPCtracker::AliTPCRow::InsertCluster(const AliTPCcluster* c, UInt_t index) {
//-----------------------------------------------------------------------
// Insert a cluster into this pad row in accordence with its y-coordinate
//-----------------------------------------------------------------------
- if (fN==kMAXCLUSTER) {
+ if (fN==kMaxClusterPerRow) {
cerr<<"AliTPCRow::InsertCluster(): Too many clusters !\n"; return;
}
if (fN==0) {fIndex[0]=index; fClusters[fN++]=c; return;}
SetdEdx(dedx);
}
-//_____________________________________________________________________________
-void AliTPCtracker::AliTPCseed::UseClusters(AliTPCClustersArray *ca, Int_t n) {
- //-----------------------------------------------------------------
- // This function marks clusters associated with this track.
- //-----------------------------------------------------------------
- Int_t nc=GetNumberOfClusters();
- Int_t sec,row,ncl;
-
- for (Int_t i=n; i<nc; i++) {
- GetCluster(i,sec,row,ncl);
- AliTPCClustersRow *clrow=ca->GetRow(sec,row);
- AliTPCcluster *c=(AliTPCcluster*)(*clrow)[ncl];
- c->Use();
- }
-}
//
// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
//-------------------------------------------------------
-
+#include "AliTracker.h"
#include "AliTPCtrack.h"
-#include "AliTPCParam.h"
+#include "AliTPCClustersArray.h"
-#define kMAXCLUSTER 2500
+#include "AliTPCreco.h"
class TFile;
+class AliTPCParam;
-class AliTPCtracker {
+class AliTPCtracker : public AliTracker {
public:
- static Int_t Clusters2Tracks(const AliTPCParam *par, TFile *of);
+ AliTPCtracker():AliTracker(),fkNIS(0),fkNOS(0) {}
+ AliTPCtracker(const AliTPCParam *par);
+ ~AliTPCtracker();
+
+ Int_t ReadSeeds(const TFile *in);
+
+ void LoadInnerSectors();
+ void UnloadInnerSectors();
+ void LoadOuterSectors();
+ void UnloadOuterSectors();
+ AliCluster *GetCluster(Int_t index) const;
+ Int_t Clusters2Tracks(const TFile *in, TFile *out);
+ Int_t PropagateBack(const TFile *in, TFile *out);
+
+ virtual void CookLabel(AliKalmanTrack *t,Float_t wrong) const;
+
+public:
//**************** Internal tracker class **********************
class AliTPCRow {
public:
const AliTPCcluster* operator[](Int_t i) const {return fClusters[i];}
UInt_t GetIndex(Int_t i) const {return fIndex[i];}
Int_t Find(Double_t y) const;
+ void SetX(Double_t x) {fX=x;}
+ Double_t GetX() const {return fX;}
private:
- unsigned fN; //number of clusters
- const AliTPCcluster *fClusters[kMAXCLUSTER]; //pointers to clusters
- UInt_t fIndex[kMAXCLUSTER]; //indeces of clusters
+ Int_t fN; //number of clusters
+ const AliTPCcluster *fClusters[kMaxClusterPerRow]; //pointers to clusters
+ UInt_t fIndex[kMaxClusterPerRow]; //indeces of clusters
+ Double_t fX; //X-coordinate of this row
private:
AliTPCRow(const AliTPCRow& r); //dummy copy constructor
class AliTPCSector {
public:
AliTPCSector() { fN=0; fRow = 0; }
- virtual ~AliTPCSector() { delete[] fRow; }
- static void SetParam(const AliTPCParam *p) { fgParam=p; }
+ ~AliTPCSector() { delete[] fRow; }
AliTPCRow& operator[](Int_t i) const { return *(fRow+i); }
Int_t GetNRows() const { return fN; }
- virtual Double_t GetX(Int_t l) const = 0;
- virtual Double_t GetMaxY(Int_t l) const = 0;
- virtual Double_t GetAlpha() const = 0;
- virtual Double_t GetAlphaShift() const = 0;
- virtual Int_t GetRowNumber(Double_t x) const = 0;
- virtual Double_t GetPadPitchWidth() const = 0;
-
- protected:
- unsigned fN; //number of pad rows
+ void Setup(const AliTPCParam *par, Int_t flag);
+ Double_t GetX(Int_t l) const {return fRow[l].GetX();}
+ Double_t GetMaxY(Int_t l) const {
+ return GetX(l)*TMath::Tan(0.5*GetAlpha());
+ }
+ Double_t GetAlpha() const {return fAlpha;}
+ Double_t GetAlphaShift() const {return fAlphaShift;}
+ Int_t GetRowNumber(Double_t x) const {
+ //return pad row number for this x
+ Double_t r=fRow[fN-1].GetX();
+ if (x > r) return fN;
+ r=fRow[0].GetX();
+ if (x < r) return -1;
+ return Int_t((x-r)/fPadPitchLength + 0.5);
+ }
+ Double_t GetPadPitchWidth() const {return fPadPitchWidth;}
+ Double_t GetPadPitchLength() const {return fPadPitchLength;}
+
+ private:
+ Int_t fN; //number of pad rows
AliTPCRow *fRow; //array of pad rows
- static const AliTPCParam *fgParam; //TPC parameters
+ Double_t fAlpha; //opening angle
+ Double_t fAlphaShift; //shift angle;
+ Double_t fPadPitchWidth; //pad pitch width
+ Double_t fPadPitchLength; //pad pitch length
private:
AliTPCSector(const AliTPCSector &s); //dummy copy contructor
AliTPCSector& operator=(const AliTPCSector &s);//dummy assignment operator
};
-//**************** Internal tracker class **********************
- class AliTPCSSector : public AliTPCSector {
- public:
- AliTPCSSector();
- Double_t GetX(Int_t l) const { return fgParam->GetPadRowRadiiLow(l); }
- Double_t GetMaxY(Int_t l) const {
- return GetX(l)*TMath::Tan(0.5*GetAlpha());
- }
- Double_t GetAlpha() const {return fgParam->GetInnerAngle();}
- Double_t GetAlphaShift() const {return fgParam->GetInnerAngleShift();}
- Double_t GetPadPitchWidth() const {
- return fgParam->GetInnerPadPitchWidth();
- }
- Int_t GetRowNumber(Double_t x) const;
- };
-
-//**************** Internal tracker class **********************
- class AliTPCLSector : public AliTPCSector {
- public:
- AliTPCLSector();
- Double_t GetX(Int_t l) const { return fgParam->GetPadRowRadiiUp(l); }
- Double_t GetMaxY(Int_t l) const {
- return GetX(l)*TMath::Tan(0.5*GetAlpha());
- }
- Double_t GetAlpha() const {return fgParam->GetOuterAngle();}
- Double_t GetAlphaShift() const {return fgParam->GetOuterAngleShift();}
- Double_t GetPadPitchWidth() const {
- return fgParam->GetOuterPadPitchWidth();
- }
- Int_t GetRowNumber(Double_t x) const;
- };
-
//**************** Internal tracker class **********************
class AliTPCseed : public AliTPCtrack {
public:
AliTPCseed():AliTPCtrack(){}
+ AliTPCseed(const AliTPCtrack &t):AliTPCtrack(t){}
+ AliTPCseed(const AliKalmanTrack &t, Double_t a):AliTPCtrack(t,a){}
AliTPCseed(UInt_t index, const Double_t xx[5],
const Double_t cc[15], Double_t xr, Double_t alpha):
AliTPCtrack(index, xx, cc, xr, alpha) {}
q *= TMath::Sqrt((1-s*s)/(1+t*t));
fdEdxSample[i]=q;
}
- void UseClusters(AliTPCClustersArray *ca, Int_t n=0);
void CookdEdx(Double_t low=0.05, Double_t up=0.70);
private:
};
private:
- static Int_t
- FindProlongation(AliTPCseed& t,const AliTPCSector *sec,Int_t s,Int_t rf=0);
- static void
- MakeSeeds(TObjArray &sd,const AliTPCSector *s,Int_t max,Int_t i1,Int_t i2);
-};
+ void MakeSeeds(Int_t i1, Int_t i2);
+ Int_t FollowProlongation(AliTPCseed& t, Int_t rf=0);
+ Int_t FollowBackProlongation(AliTPCseed &s, const AliTPCtrack &t);
+
+ AliTPCtracker(const AliTPCtracker& r); //dummy copy constructor
+ AliTPCtracker &operator=(const AliTPCtracker& r);//dummy assignment operator
+
+ const Int_t fkNIS; //number of inner sectors
+ AliTPCSector *fInnerSec; //array of inner sectors;
+ const Int_t fkNOS; //number of outer sectors
+ AliTPCSector *fOuterSec; //array of outer sectors;
-inline AliTPCtracker::AliTPCSSector::AliTPCSSector() {
- //default constructor
- if (!fgParam) {
- fprintf(stderr,"AliTPCSSector: parameters are not set !\n");
- return;
- }
- fN=fgParam->GetNRowLow();
- fRow=new AliTPCRow[fN];
-}
-
-inline
-Int_t AliTPCtracker::AliTPCSSector::GetRowNumber(Double_t x) const {
- //return pad row number for this x
- Double_t r=fgParam->GetPadRowRadiiLow(fgParam->GetNRowLow()-1);
- if (x > r) return fgParam->GetNRowLow();
- r=fgParam->GetPadRowRadiiLow(0);
- if (x < r) return -1;
- return Int_t((x-r)/fgParam->GetInnerPadPitchLength() + 0.5);
-}
-
-inline AliTPCtracker::AliTPCLSector::AliTPCLSector(){
- //default constructor
- if (!fgParam) {
- fprintf(stderr,"AliTPCLSector: parameters are not set !\n");
- return;
- }
- fN=fgParam->GetNRowUp();
- fRow=new AliTPCRow[fN];
-}
-
-inline
-Int_t AliTPCtracker::AliTPCLSector::GetRowNumber(Double_t x) const {
- //return pad row number for this x
- Double_t r=fgParam->GetPadRowRadiiUp(fgParam->GetNRowUp()-1);
- if (x > r) return fgParam->GetNRowUp();
- r=fgParam->GetPadRowRadiiUp(0);
- if (x < r) return -1;
- return Int_t((x-r)/fgParam->GetOuterPadPitchLength() + 0.5);
-}
-
-//-----------------------------------------------------------------
+ Int_t fN; //number of loaded sectors
+ AliTPCSector *fSectors; //pointer to loaded sectors;
+
+ AliTPCClustersArray fClustersArray; //array of TPC clusters
+ TObjArray *fSeeds; //array of track seeds
+};
#endif
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff