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
13 //-----------------------------------------------------------------------------
15 //----- Root headers ---------
17 //---- AliRoot headers -------
19 #include "AliGausCorr.h"
21 #include "AliTPCkineGrid.h"
22 #include "AliTPCtrack.h"
23 #include "AliTrackReference.h"
24 //----------------------------
26 class AliTPCtrackerParam {
27 /////////////////////////////////////////////////////////////////////////
29 // This class builds AliTPCtrack objects from generated tracks to feed
30 // ITS tracking (V2). The AliTPCtrack is built from its first hit in
31 // the TPC. The track is assigned a Kalman-like covariance matrix
32 // depending on its pT and pseudorapidity and track parameters are
33 // smeared according to this covariance matrix.
34 // Output file contains sorted tracks, ready for matching with ITS.
36 // See implementation file for more details.
39 // Origin: Andrea Dainese, Padova - e-mail: andrea.dainese@pd.infn.it
41 /////////////////////////////////////////////////////////////////////////
43 AliTPCtrackerParam(const Int_t coll=0,const Double_t Bz=0.4,const Int_t n=1);
44 virtual ~AliTPCtrackerParam();
46 // this function performs the parameterized tracking
47 Int_t BuildTPCtracks(const TFile *inp, TFile *out);
49 // these functions are used to create a DB of cov. matrices,
50 // including regularization, efficiencies and dE/dx
51 void AllGeantTracks() { fSelAndSmear=kFALSE; return; }
52 void AnalyzedEdx(const Char_t *outName,Int_t pdg);
53 void AnalyzePulls(const Char_t *outName);
54 void AnalyzeResolutions(Int_t pdg);
55 void CompareTPCtracks(const Char_t *galiceName="galice.root",
56 const Char_t *trkGeaName="AliTPCtracksGeant.root",
57 const Char_t *trkKalName="AliTPCtracksSorted.root",
58 const Char_t *covmatName="CovMatrix.root",
59 const Char_t *tpceffasciiName="TPCeff.dat",
60 const Char_t *tpceffrootName="TPCeff.root");
61 void DrawEffs(const Char_t *inName,Int_t pdg=211);
62 void DrawPulls(const Char_t *inName,Int_t pdg=211,Int_t par=0);
64 void MergeEvents(Int_t evFirst=1,Int_t evLast=1);
65 void RegularizeCovMatrix(const Char_t *outName,Int_t pdg);
68 //********* Internal class definition *******
69 class AliTPCtrackParam : public AliTPCtrack {
71 AliTPCtrackParam():AliTPCtrack(){}
72 AliTPCtrackParam(const AliTPCtrack &t):AliTPCtrack(t){}
74 void AssignMass(Double_t mass) {SetMass(mass); return;}
79 //********* end of internal class ***********
81 //********* Internal class definition *******
82 class AliTPCseedGeant : public TObject {
84 AliTPCseedGeant(Double_t x,Double_t y,Double_t z,
85 Double_t px,Double_t py,Double_t pz,
94 Double_t a = TMath::ATan2(y,x)*180./TMath::Pi();
96 fSector = (Int_t)(a/20.);
97 fAlpha = 10.+20.*fSector;
99 fAlpha *= TMath::Pi();
101 Int_t GetLabel() { return fLabel; }
102 Double_t GetAlpha() { return fAlpha; }
103 Double_t GetXL() { return fXg*TMath::Cos(fAlpha)+fYg*TMath::Sin(fAlpha); }
104 Double_t GetYL() { return -fXg*TMath::Sin(fAlpha)+fYg*TMath::Cos(fAlpha); }
105 Double_t GetZL() { return fZg; }
106 Double_t GetPx() { return fPx; }
107 Double_t GetPy() { return fPy; }
108 Double_t GetPz() { return fPz; }
109 Double_t GetPt() { return TMath::Sqrt(fPx*fPx+fPy*fPy); }
110 Double_t GetEta() { return -TMath::Log(TMath::Tan(0.25*TMath::Pi()-0.5*TMath::ATan(fPz/GetPt()))); }
111 void SetLabel(Int_t lab) { fLabel=lab; return; }
112 Bool_t InTPCAcceptance() {
113 if(TMath::Abs(GetZL()+(244.-GetXL())*fPz/GetPt())>252.) return kFALSE;
128 //******* end of internal class ****************
131 Int_t fNevents; // number of events in the file to be processed
132 Double_t fBz; // value of the z component of L3 field (Tesla)
133 Int_t fColl; // collision code (0: PbPb6000; 1: pp)
134 Bool_t fSelAndSmear; // if kFALSE returns GEANT tracks
136 TString fDBfileName; // DataBase file name
138 AliTPCtrack fTrack; // current track
140 TTree *fCovTree; // tree with regularized cov matrices
141 // for the current track
143 AliTPCkineGrid *fDBgrid; // grid for the cov matrix look-up table
144 AliTPCkineGrid fDBgridPi; // " for pions
145 AliTPCkineGrid fDBgridKa; // " for kaons
146 AliTPCkineGrid fDBgridPr; // " for protons
147 AliTPCkineGrid fDBgridEl; // " for electrons
148 AliTPCkineGrid fDBgridMu; // " for muons
150 AliTPCkineGrid *fEff; // TPC efficiencies for the current track
151 AliTPCkineGrid fEffPi; // " pions
152 AliTPCkineGrid fEffKa; // " kaons
153 AliTPCkineGrid fEffPr; // " protons
154 AliTPCkineGrid fEffEl; // " electrons
155 AliTPCkineGrid fEffMu; // " muons
157 AliTPCkineGrid *fPulls; // pulls for the current track
158 AliTPCkineGrid fPullsPi[5]; // " pions
159 AliTPCkineGrid fPullsKa[5]; // " muons
160 AliTPCkineGrid fPullsPr[5]; // " protons
161 AliTPCkineGrid fPullsEl[5]; // " electrons
162 AliTPCkineGrid fPullsMu[5]; // " muons
164 TMatrixD *fRegPar; // regularization parameters for the curr. track
165 TMatrixD fRegParPi; // " for pions
166 TMatrixD fRegParKa; // " for kaons
167 TMatrixD fRegParPr; // " for protons
168 TMatrixD fRegParEl; // " for electrons
169 TMatrixD fRegParMu; // " for muons
171 AliTPCkineGrid *fdEdxMean; // dEdx mean for the current track
172 AliTPCkineGrid fdEdxMeanPi; // " pions
173 AliTPCkineGrid fdEdxMeanKa; // " kaons
174 AliTPCkineGrid fdEdxMeanPr; // " protons
175 AliTPCkineGrid fdEdxMeanEl; // " electrons
176 AliTPCkineGrid fdEdxMeanMu; // " muons
178 AliTPCkineGrid *fdEdxRMS; // dEdx RMS for the current track
179 AliTPCkineGrid fdEdxRMSPi; // " pions
180 AliTPCkineGrid fdEdxRMSKa; // " kaons
181 AliTPCkineGrid fdEdxRMSPr; // " protons
182 AliTPCkineGrid fdEdxRMSEl; // " electrons
183 AliTPCkineGrid fdEdxRMSMu; // " muons
186 void BuildTrack(AliTPCseedGeant *s,Int_t ch);
187 Int_t CheckLabel(AliTPCseedGeant *s,Int_t nPart,
188 Double_t *ptkine,Double_t *pzkine) const;
189 void CookdEdx(Double_t pt,Double_t eta);
190 void CookTrack(Double_t pt,Double_t eta);
191 TMatrixD GetSmearingMatrix(Double_t *cc, Double_t pt,Double_t eta) const;
192 void InitializeKineGrid(Option_t *which);
193 void MakeSeedsFromHits(AliTPC *TPC,TTree *TH,TObjArray &seedArray) const;
194 void MakeSeedsFromRefs(TTree *TTR,
195 TObjArray &seedArray) const;
196 Int_t ReadAllData(const Char_t *inName);
197 Int_t ReadDBgrid(const Char_t *inName);
198 Int_t ReaddEdx(const Char_t *inName,Int_t pdg);
199 Int_t ReadEffs(const Char_t *inName);
200 Int_t ReadPulls(const Char_t *inName);
201 Int_t ReadRegParams(const Char_t *inName,Int_t pdg);
202 Bool_t SelectedTrack(Double_t pt, Double_t eta) const;
203 void SetParticle(Int_t pdg);
204 void SmearTrack(Double_t *xx,Double_t *xxsm,TMatrixD cov) const;
205 Int_t WritedEdx(const Char_t *outName,Int_t pdg);
206 Int_t WriteEffs(const Char_t *outName);
207 Int_t WritePulls(const Char_t *outName);
208 Int_t WriteRegParams(const Char_t *outName,Int_t pdg);
211 ClassDef(AliTPCtrackerParam,1) // TPC tracking parameterization class