1 #ifndef ALITPCTRACKERPARAM_H
2 #define ALITPCTRACKERPARAM_H
3 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. */
4 /* See cxx source for full Copyright notice */
9 //-----------------------------------------------------------------------------
10 // TPC Tracking Parameterization Class
12 // Origin: Andrea Dainese, Padova - e-mail: andrea.dainese@pd.infn.it
16 //-----------------------------------------------------------------------------
18 //----- Root headers ---------
20 //---- AliRoot headers -------
23 //#include "AliGausCorr.h"
24 //#include "AliMagF.h"
25 #include "AliTPCkineGrid.h"
26 #include "AliTPCtrack.h"
27 //#include "AliTrackReference.h"
28 //----------------------------
30 class AliTPCtrackerParam:
34 /////////////////////////////////////////////////////////////////////////
36 // This class builds AliTPCtrack objects from generated tracks to feed
37 // ITS tracking (V2). The AliTPCtrack is built from its first hit in
38 // the TPC. The track is assigned a Kalman-like covariance matrix
39 // depending on its pT and pseudorapidity and track parameters are
40 // smeared according to this covariance matrix.
41 // Output file contains sorted tracks, ready for matching with ITS.
43 // See implementation file for more details.
46 // Origin: Andrea Dainese, Padova - e-mail: andrea.dainese@pd.infn.it
48 /////////////////////////////////////////////////////////////////////////
50 AliTPCtrackerParam(Int_t coll=0, Double_t Bz=0.4, Int_t n=1,
51 const char* evfoldname = AliConfig::GetDefaultEventFolderName());
52 virtual ~AliTPCtrackerParam();
54 // this function performs the parameterized tracking
56 AliTPCtrackerParam(const AliTPCtrackerParam& p);
58 Int_t BuildTPCtracks(const TFile *inp, TFile *out);
59 // these functions are used to create a DB of cov. matrices,
60 // including regularization, efficiencies and dE/dx
61 void AllGeantTracks() { fSelAndSmear=kFALSE; return; }
62 void AnalyzedEdx(const Char_t *outName,Int_t pdg);
63 void AnalyzePulls(const Char_t *outName);
64 void AnalyzeResolutions(Int_t pdg);
65 void CompareTPCtracks(const Char_t *galiceName="galice.root",
66 const Char_t *trkGeaName="AliTPCtracksGeant.root",
67 const Char_t *trkKalName="AliTPCtracksSorted.root",
68 const Char_t *covmatName="CovMatrix.root",
69 const Char_t *tpceffasciiName="TPCeff.dat",
70 const Char_t *tpceffrootName="TPCeff.root");
71 void DrawEffs(const Char_t *inName,Int_t pdg=211);
72 void DrawPulls(const Char_t *inName,Int_t pdg=211,Int_t par=0);
74 void MergeEvents(Int_t evFirst=1,Int_t evLast=1);
75 void RegularizeCovMatrix(const Char_t *outName,Int_t pdg);
78 //********* Internal class definition *******
79 class AliTPCtrackParam : public AliTPCtrack {
81 AliTPCtrackParam():AliTPCtrack(){}
82 AliTPCtrackParam(const AliTPCtrack &t):AliTPCtrack(t){}
84 void AssignMass(Double_t mass) {SetMass(mass); return;}
89 //********* end of internal class ***********
91 //********* Internal class definition *******
92 class AliTPCseedGeant : public TObject {
94 AliTPCseedGeant(Double_t x=0.,Double_t y=0.,Double_t z=0.,
95 Double_t px=0.,Double_t py=0.,Double_t pz=0.,
97 Int_t GetLabel() const { return fLabel; }
98 Double_t GetAlpha() const { return fAlpha; }
99 Double_t GetXL() const
100 { return fXg*TMath::Cos(fAlpha)+fYg*TMath::Sin(fAlpha); }
101 Double_t GetYL() const
102 { return -fXg*TMath::Sin(fAlpha)+fYg*TMath::Cos(fAlpha); }
103 Double_t GetZL() const { return fZg; }
104 Double_t GetPx() const { return fPx; }
105 Double_t GetPy() const { return fPy; }
106 Double_t GetPz() const { return fPz; }
107 Double_t GetPt() const { return TMath::Sqrt(fPx*fPx+fPy*fPy); }
108 Double_t GetEta() const
109 { return -TMath::Log(TMath::Tan(0.25*TMath::Pi()-0.5*TMath::ATan(fPz/GetPt()))); }
110 void SetLabel(Int_t lab) { fLabel=lab; return; }
111 Bool_t InTPCAcceptance() const {
112 if(TMath::Abs(GetZL()+(244.-GetXL())*fPz/GetPt())>252.) return kFALSE;
117 Double_t fXg; // global x of seed
118 Double_t fYg; // global y of seed
119 Double_t fZg; // global z of seed
120 Double_t fPx; // global px
121 Double_t fPy; // global py
122 Double_t fPz; // global pz
123 Double_t fAlpha; // alpha angle
124 Int_t fLabel; // track label
125 Int_t fSector; // TPC sector
127 //******* end of internal class ****************
130 TString fEvFolderName;//! name of data folder
132 Int_t fNevents; // number of events in the file to be processed
133 Double_t fBz; // value of the z component of L3 field (Tesla)
134 Int_t fColl; // collision code (0: PbPb6000; 1: pp)
135 Bool_t fSelAndSmear; // if kFALSE returns GEANT tracks
137 TString fDBfileName; // DataBase file name
139 AliTPCtrack fTrack; // current track
141 TTree *fCovTree; // tree with regularized cov matrices
142 // for the current track
144 AliTPCkineGrid *fDBgrid; // grid for the cov matrix look-up table
145 AliTPCkineGrid fDBgridPi; // " for pions
146 AliTPCkineGrid fDBgridKa; // " for kaons
147 AliTPCkineGrid fDBgridPr; // " for protons
148 AliTPCkineGrid fDBgridEl; // " for electrons
149 AliTPCkineGrid fDBgridMu; // " for muons
151 AliTPCkineGrid *fEff; // TPC efficiencies for the current track
152 AliTPCkineGrid fEffPi; // " pions
153 AliTPCkineGrid fEffKa; // " kaons
154 AliTPCkineGrid fEffPr; // " protons
155 AliTPCkineGrid fEffEl; // " electrons
156 AliTPCkineGrid fEffMu; // " muons
158 AliTPCkineGrid *fPulls; // pulls for the current track
159 AliTPCkineGrid fPullsPi[5]; // " pions
160 AliTPCkineGrid fPullsKa[5]; // " muons
161 AliTPCkineGrid fPullsPr[5]; // " protons
162 AliTPCkineGrid fPullsEl[5]; // " electrons
163 AliTPCkineGrid fPullsMu[5]; // " muons
165 TMatrixD *fRegPar; // regularization parameters for the curr. track
166 TMatrixD fRegParPi; // " for pions
167 TMatrixD fRegParKa; // " for kaons
168 TMatrixD fRegParPr; // " for protons
169 TMatrixD fRegParEl; // " for electrons
170 TMatrixD fRegParMu; // " for muons
172 AliTPCkineGrid *fdEdxMean; // dEdx mean for the current track
173 AliTPCkineGrid fdEdxMeanPi; // " pions
174 AliTPCkineGrid fdEdxMeanKa; // " kaons
175 AliTPCkineGrid fdEdxMeanPr; // " protons
176 AliTPCkineGrid fdEdxMeanEl; // " electrons
177 AliTPCkineGrid fdEdxMeanMu; // " muons
179 AliTPCkineGrid *fdEdxRMS; // dEdx RMS for the current track
180 AliTPCkineGrid fdEdxRMSPi; // " pions
181 AliTPCkineGrid fdEdxRMSKa; // " kaons
182 AliTPCkineGrid fdEdxRMSPr; // " protons
183 AliTPCkineGrid fdEdxRMSEl; // " electrons
184 AliTPCkineGrid fdEdxRMSMu; // " muons
187 void BuildTrack(AliTPCseedGeant *s,Int_t ch);
188 Int_t CheckLabel(AliTPCseedGeant *s,Int_t nPart,
189 Double_t *ptkine,Double_t *pzkine) const;
190 void CookdEdx(Double_t pt,Double_t eta);
191 void CookTrack(Double_t pt,Double_t eta);
192 TMatrixD GetSmearingMatrix(Double_t *cc, Double_t pt,Double_t eta) const;
193 void InitializeKineGrid(Option_t *which);
194 void MakeSeedsFromHits(AliTPC *tpc,TTree *th,TObjArray &seedArray) const;
195 void MakeSeedsFromRefs(TTree *ttr,
196 TObjArray &seedArray) const;
197 Int_t ReadAllData(const Char_t *inName);
198 Int_t ReadDBgrid(const Char_t *inName);
199 Int_t ReaddEdx(const Char_t *inName,Int_t pdg);
200 Int_t ReadEffs(const Char_t *inName);
201 Int_t ReadPulls(const Char_t *inName);
202 Int_t ReadRegParams(const Char_t *inName,Int_t pdg);
203 Bool_t SelectedTrack(Double_t pt, Double_t eta) const;
204 void SetParticle(Int_t pdg);
205 void SmearTrack(Double_t *xx,Double_t *xxsm,TMatrixD cov) const;
206 Int_t WritedEdx(const Char_t *outName,Int_t pdg);
207 Int_t WriteEffs(const Char_t *outName);
208 Int_t WritePulls(const Char_t *outName);
209 Int_t WriteRegParams(const Char_t *outName,Int_t pdg);
212 ClassDef(AliTPCtrackerParam,1) // TPC tracking parameterization class