]>
Commit | Line | Data |
---|---|---|
46d29e70 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | /* | |
17 | $Log$ | |
79e94bf8 | 18 | Revision 1.26 2003/04/10 10:36:54 hristov |
19 | Code for unified TPC/TRD tracking (S.Radomski) | |
20 | ||
9c9d2487 | 21 | Revision 1.25 2003/03/19 17:14:11 hristov |
22 | Load/UnloadClusters added to the base class and the derived classes changed correspondingly. Possibility to give 2 input files for ITS and TPC tracks in PropagateBack. TRD tracker uses fEventN from the base class (T.Kuhr) | |
23 | ||
e24ea474 | 24 | Revision 1.24 2003/02/19 09:02:28 hristov |
25 | Track time measurement (S.Radomski) | |
26 | ||
0d5b5c27 | 27 | Revision 1.23 2003/02/10 14:06:10 cblume |
28 | Add tracking without tilted pads as option | |
29 | ||
b8dc2353 | 30 | Revision 1.22 2003/01/30 15:19:58 cblume |
31 | New set of parameters | |
32 | ||
a9814c09 | 33 | Revision 1.21 2003/01/27 16:34:49 cblume |
34 | Update of tracking by Sergei and Chuncheng | |
35 | ||
b3a5a838 | 36 | Revision 1.20 2002/11/07 15:52:09 cblume |
37 | Update of tracking code for tilted pads | |
38 | ||
fd621f36 | 39 | Revision 1.19 2002/10/22 15:53:08 alibrary |
40 | Introducing Riostream.h | |
41 | ||
a2cb5b3d | 42 | Revision 1.18 2002/10/14 14:57:44 hristov |
43 | Merging the VirtualMC branch to the main development branch (HEAD) | |
44 | ||
b9d0a01d | 45 | Revision 1.14.6.2 2002/07/24 10:09:31 alibrary |
46 | Updating VirtualMC | |
47 | ||
48 | Revision 1.17 2002/06/13 12:09:58 hristov | |
49 | Minor corrections | |
50 | ||
ca6c93d6 | 51 | Revision 1.16 2002/06/12 09:54:36 cblume |
52 | Update of tracking code provided by Sergei | |
53 | ||
0a29d0f1 | 54 | Revision 1.14 2001/11/14 10:50:46 cblume |
55 | Changes in digits IO. Add merging of summable digits | |
56 | ||
abaf1f1d | 57 | Revision 1.13 2001/05/30 12:17:47 hristov |
58 | Loop variables declared once | |
59 | ||
3ab6f951 | 60 | Revision 1.12 2001/05/28 17:07:58 hristov |
61 | Last minute changes; ExB correction in AliTRDclusterizerV1; taking into account of material in G10 TEC frames and material between TEC planes (C.Blume,S.Sedykh) | |
62 | ||
71d9fa7b | 63 | Revision 1.8 2000/12/20 13:00:44 cblume |
64 | Modifications for the HP-compiler | |
65 | ||
d1b06c24 | 66 | Revision 1.7 2000/12/08 16:07:02 cblume |
67 | Update of the tracking by Sergei | |
68 | ||
bbf92647 | 69 | Revision 1.6 2000/11/30 17:38:08 cblume |
70 | Changes to get in line with new STEER and EVGEN | |
71 | ||
1819f4bb | 72 | Revision 1.5 2000/11/14 14:40:27 cblume |
73 | Correction for the Sun compiler (kTRUE and kFALSE) | |
74 | ||
57527628 | 75 | Revision 1.4 2000/11/10 14:57:52 cblume |
76 | Changes in the geometry constants for the DEC compiler | |
77 | ||
dd56b762 | 78 | Revision 1.3 2000/10/15 23:40:01 cblume |
79 | Remove AliTRDconst | |
80 | ||
0e9c2ad5 | 81 | Revision 1.2 2000/10/06 16:49:46 cblume |
82 | Made Getters const | |
83 | ||
46d29e70 | 84 | Revision 1.1.2.2 2000/10/04 16:34:58 cblume |
85 | Replace include files by forward declarations | |
86 | ||
87 | Revision 1.1.2.1 2000/09/22 14:47:52 cblume | |
88 | Add the tracking code | |
89 | ||
90 | */ | |
bbf92647 | 91 | |
a2cb5b3d | 92 | #include <Riostream.h> |
46d29e70 | 93 | |
94 | #include <TFile.h> | |
46d29e70 | 95 | #include <TBranch.h> |
5443e65e | 96 | #include <TTree.h> |
97 | #include <TObjArray.h> | |
46d29e70 | 98 | |
46d29e70 | 99 | #include "AliTRDgeometry.h" |
5443e65e | 100 | #include "AliTRDparameter.h" |
101 | #include "AliTRDgeometryDetail.h" | |
46d29e70 | 102 | #include "AliTRDcluster.h" |
103 | #include "AliTRDtrack.h" | |
79e94bf8 | 104 | #include "AliTRDPartID.h" |
5443e65e | 105 | #include "../TPC/AliTPCtrack.h" |
46d29e70 | 106 | |
107 | #include "AliTRDtracker.h" | |
108 | ||
109 | ClassImp(AliTRDtracker) | |
110 | ||
5443e65e | 111 | const Float_t AliTRDtracker::fSeedDepth = 0.5; |
112 | const Float_t AliTRDtracker::fSeedStep = 0.10; | |
113 | const Float_t AliTRDtracker::fSeedGap = 0.25; | |
114 | ||
115 | const Float_t AliTRDtracker::fMaxSeedDeltaZ12 = 40.; | |
116 | const Float_t AliTRDtracker::fMaxSeedDeltaZ = 25.; | |
117 | const Float_t AliTRDtracker::fMaxSeedC = 0.0052; | |
118 | const Float_t AliTRDtracker::fMaxSeedTan = 1.2; | |
119 | const Float_t AliTRDtracker::fMaxSeedVertexZ = 150.; | |
a819a5f7 | 120 | |
5443e65e | 121 | const Double_t AliTRDtracker::fSeedErrorSY = 0.2; |
122 | const Double_t AliTRDtracker::fSeedErrorSY3 = 2.5; | |
123 | const Double_t AliTRDtracker::fSeedErrorSZ = 0.1; | |
bbf92647 | 124 | |
5443e65e | 125 | const Float_t AliTRDtracker::fMinClustersInSeed = 0.7; |
bbf92647 | 126 | |
5443e65e | 127 | const Float_t AliTRDtracker::fMinClustersInTrack = 0.5; |
128 | const Float_t AliTRDtracker::fMinFractionOfFoundClusters = 0.8; | |
bbf92647 | 129 | |
a9814c09 | 130 | const Float_t AliTRDtracker::fSkipDepth = 0.3; |
5443e65e | 131 | const Float_t AliTRDtracker::fLabelFraction = 0.8; |
132 | const Float_t AliTRDtracker::fWideRoad = 20.; | |
133 | ||
b8dc2353 | 134 | const Double_t AliTRDtracker::fMaxChi2 = 12.; |
a819a5f7 | 135 | |
9c9d2487 | 136 | const Int_t AliTRDtracker::kFirstPlane = 5; |
137 | const Int_t AliTRDtracker::kLastPlane = 17; | |
138 | ||
bbf92647 | 139 | |
46d29e70 | 140 | //____________________________________________________________________ |
9c9d2487 | 141 | AliTRDtracker::AliTRDtracker(const TFile *geomfile):AliTracker() |
46d29e70 | 142 | { |
5443e65e | 143 | // |
144 | // Main constructor | |
145 | // | |
46d29e70 | 146 | |
5443e65e | 147 | Float_t fTzero = 0; |
b8dc2353 | 148 | |
5443e65e | 149 | fAddTRDseeds = kFALSE; |
5443e65e | 150 | fGeom = NULL; |
b8dc2353 | 151 | fNoTilt = kFALSE; |
5443e65e | 152 | |
153 | TDirectory *savedir=gDirectory; | |
154 | TFile *in=(TFile*)geomfile; | |
155 | if (!in->IsOpen()) { | |
156 | printf("AliTRDtracker::AliTRDtracker(): geometry file is not open!\n"); | |
157 | printf(" DETAIL TRD geometry and DEFAULT TRD parameter will be used\n"); | |
158 | } | |
159 | else { | |
160 | in->cd(); | |
79e94bf8 | 161 | // in->ls(); |
5443e65e | 162 | fGeom = (AliTRDgeometry*) in->Get("TRDgeometry"); |
163 | fPar = (AliTRDparameter*) in->Get("TRDparameter"); | |
79e94bf8 | 164 | // fGeom->Dump(); |
5443e65e | 165 | } |
46d29e70 | 166 | |
5443e65e | 167 | if(fGeom) { |
168 | // fTzero = geo->GetT0(); | |
b8dc2353 | 169 | printf("Found geometry version %d on file \n", fGeom->IsVersion()); |
5443e65e | 170 | } |
171 | else { | |
172 | printf("AliTRDtracker::AliTRDtracker(): cann't find TRD geometry!\n"); | |
173 | printf(" DETAIL TRD geometry and DEFAULT TRD parameter will be used\n"); | |
174 | fGeom = new AliTRDgeometryDetail(); | |
175 | fPar = new AliTRDparameter(); | |
176 | } | |
177 | ||
178 | savedir->cd(); | |
46d29e70 | 179 | |
5443e65e | 180 | |
181 | // fGeom->SetT0(fTzero); | |
0a29d0f1 | 182 | |
46d29e70 | 183 | fNclusters = 0; |
184 | fClusters = new TObjArray(2000); | |
185 | fNseeds = 0; | |
5443e65e | 186 | fSeeds = new TObjArray(2000); |
46d29e70 | 187 | fNtracks = 0; |
5443e65e | 188 | fTracks = new TObjArray(1000); |
a819a5f7 | 189 | |
5443e65e | 190 | for(Int_t geom_s = 0; geom_s < kTRACKING_SECTORS; geom_s++) { |
191 | Int_t tr_s = CookSectorIndex(geom_s); | |
192 | fTrSec[tr_s] = new AliTRDtrackingSector(fGeom, geom_s, fPar); | |
193 | } | |
a819a5f7 | 194 | |
b8dc2353 | 195 | Float_t tilt_angle = TMath::Abs(fPar->GetTiltingAngle()); |
196 | if(tilt_angle < 0.1) { | |
197 | fNoTilt = kTRUE; | |
198 | } | |
199 | ||
200 | fSY2corr = 0.2; | |
201 | fSZ2corr = 120.; | |
202 | ||
203 | if(fNoTilt && (tilt_angle > 0.1)) fSY2corr = fSY2corr + tilt_angle * 0.05; | |
204 | ||
bbf92647 | 205 | |
5443e65e | 206 | // calculate max gap on track |
a819a5f7 | 207 | |
5443e65e | 208 | Double_t dxAmp = (Double_t) fGeom->CamHght(); // Amplification region |
209 | Double_t dxDrift = (Double_t) fGeom->CdrHght(); // Drift region | |
a819a5f7 | 210 | |
5443e65e | 211 | Double_t dx = (Double_t) fPar->GetTimeBinSize(); |
212 | Int_t tbAmp = fPar->GetTimeBefore(); | |
213 | Int_t maxAmp = (Int_t) ((dxAmp+0.000001)/dx); | |
b3a5a838 | 214 | if(kTRUE) maxAmp = 0; // intentional until we change the parameter class |
5443e65e | 215 | Int_t tbDrift = fPar->GetTimeMax(); |
216 | Int_t maxDrift = (Int_t) ((dxDrift+0.000001)/dx); | |
a819a5f7 | 217 | |
5443e65e | 218 | tbDrift = TMath::Min(tbDrift,maxDrift); |
219 | tbAmp = TMath::Min(tbAmp,maxAmp); | |
46d29e70 | 220 | |
5443e65e | 221 | fTimeBinsPerPlane = tbAmp + tbDrift; |
222 | fMaxGap = (Int_t) (fTimeBinsPerPlane * fGeom->Nplan() * fSkipDepth); | |
46d29e70 | 223 | |
5443e65e | 224 | fVocal = kFALSE; |
0a29d0f1 | 225 | |
9c9d2487 | 226 | |
227 | // Barrel Tracks [SR, 03.04.2003] | |
228 | ||
229 | fBarrelFile = 0; | |
230 | fBarrelTree = 0; | |
231 | fBarrelArray = 0; | |
232 | fBarrelTrack = 0; | |
233 | ||
234 | savedir->cd(); | |
5443e65e | 235 | } |
46d29e70 | 236 | |
5443e65e | 237 | //___________________________________________________________________ |
238 | AliTRDtracker::~AliTRDtracker() | |
46d29e70 | 239 | { |
5443e65e | 240 | delete fClusters; |
241 | delete fTracks; | |
242 | delete fSeeds; | |
243 | delete fGeom; | |
244 | delete fPar; | |
0a29d0f1 | 245 | |
5443e65e | 246 | for(Int_t geom_s = 0; geom_s < kTRACKING_SECTORS; geom_s++) { |
247 | delete fTrSec[geom_s]; | |
248 | } | |
249 | } | |
46d29e70 | 250 | |
9c9d2487 | 251 | //_____________________________________________________________________ |
252 | ||
253 | void AliTRDtracker::SetBarrelTree(const char *mode) { | |
254 | // | |
255 | // | |
256 | // | |
257 | ||
258 | if (!IsStoringBarrel()) return; | |
259 | ||
260 | TDirectory *sav = gDirectory; | |
261 | if (!fBarrelFile) fBarrelFile = new TFile("AliBarrelTracks.root", "UPDATE"); | |
262 | ||
263 | char buff[40]; | |
264 | sprintf(buff, "BarrelTRD_%d_%s", GetEventNumber(), mode); | |
265 | ||
266 | fBarrelFile->cd(); | |
267 | fBarrelTree = new TTree(buff, "Barrel TPC tracks"); | |
268 | ||
269 | Int_t nRefs = kLastPlane - kFirstPlane + 1; | |
270 | ||
271 | if (!fBarrelArray) fBarrelArray = new TClonesArray("AliBarrelTrack", nRefs); | |
272 | for(Int_t i=0; i<nRefs; i++) new((*fBarrelArray)[i]) AliBarrelTrack(); | |
273 | ||
274 | fBarrelTree->Branch("tracks", &fBarrelArray); | |
275 | sav->cd(); | |
276 | } | |
277 | ||
278 | //_____________________________________________________________________ | |
279 | ||
280 | void AliTRDtracker::StoreBarrelTrack(AliTRDtrack *ps, Int_t refPlane, Int_t isIn) { | |
281 | // | |
282 | // | |
283 | // | |
284 | ||
285 | if (!IsStoringBarrel()) return; | |
286 | ||
287 | static Int_t nClusters; | |
288 | static Int_t nWrong; | |
289 | static Double_t chi2; | |
290 | static Int_t index; | |
291 | static Bool_t wasLast = kTRUE; | |
292 | ||
293 | Int_t newClusters, newWrong; | |
294 | Double_t newChi2; | |
295 | ||
296 | if (wasLast) { | |
297 | ||
298 | fBarrelArray->Clear(); | |
299 | nClusters = nWrong = 0; | |
300 | chi2 = 0.0; | |
301 | index = 0; | |
302 | wasLast = kFALSE; | |
303 | } | |
304 | ||
305 | fBarrelTrack = (AliBarrelTrack*)(*fBarrelArray)[index++]; | |
306 | ps->GetBarrelTrack(fBarrelTrack); | |
307 | ||
308 | newClusters = ps->GetNumberOfClusters() - nClusters; | |
309 | newWrong = ps->GetNWrong() - nWrong; | |
310 | newChi2 = ps->GetChi2() - chi2; | |
311 | ||
312 | nClusters = ps->GetNumberOfClusters(); | |
313 | nWrong = ps->GetNWrong(); | |
314 | chi2 = ps->GetChi2(); | |
315 | ||
316 | if (refPlane != kLastPlane) { | |
317 | fBarrelTrack->SetNClusters(newClusters, newChi2); | |
318 | fBarrelTrack->SetNWrongClusters(newWrong); | |
319 | } else { | |
320 | wasLast = kTRUE; | |
321 | } | |
322 | ||
323 | fBarrelTrack->SetRefPlane(refPlane, isIn); | |
324 | } | |
325 | ||
326 | //_____________________________________________________________________ | |
327 | ||
328 | Bool_t AliTRDtracker::AdjustSector(AliTRDtrack *track) { | |
329 | // | |
330 | // Rotates the track when necessary | |
331 | // | |
332 | ||
333 | Double_t alpha = AliTRDgeometry::GetAlpha(); | |
334 | Double_t y = track->GetY(); | |
335 | Double_t ymax = track->GetX()*TMath::Tan(0.5*alpha); | |
336 | ||
337 | Int_t ns = AliTRDgeometry::kNsect; | |
338 | //Int_t s=Int_t(track->GetAlpha()/alpha)%ns; | |
339 | ||
340 | if (y > ymax) { | |
341 | //s = (s+1) % ns; | |
342 | if (!track->Rotate(alpha)) return kFALSE; | |
343 | } else if (y <-ymax) { | |
344 | //s = (s-1+ns) % ns; | |
345 | if (!track->Rotate(-alpha)) return kFALSE; | |
346 | } | |
347 | ||
348 | return kTRUE; | |
349 | } | |
350 | ||
46d29e70 | 351 | //_____________________________________________________________________ |
5443e65e | 352 | inline Double_t f1trd(Double_t x1,Double_t y1, |
a9814c09 | 353 | Double_t x2,Double_t y2, |
354 | Double_t x3,Double_t y3) | |
46d29e70 | 355 | { |
0a29d0f1 | 356 | // |
46d29e70 | 357 | // Initial approximation of the track curvature |
0a29d0f1 | 358 | // |
46d29e70 | 359 | Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1); |
360 | Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)- | |
361 | (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2)); | |
362 | Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)- | |
363 | (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1)); | |
364 | ||
365 | Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b); | |
366 | ||
367 | return -xr*yr/sqrt(xr*xr+yr*yr); | |
368 | } | |
369 | ||
370 | //_____________________________________________________________________ | |
5443e65e | 371 | inline Double_t f2trd(Double_t x1,Double_t y1, |
a9814c09 | 372 | Double_t x2,Double_t y2, |
373 | Double_t x3,Double_t y3) | |
46d29e70 | 374 | { |
0a29d0f1 | 375 | // |
5443e65e | 376 | // Initial approximation of the track curvature times X coordinate |
377 | // of the center of curvature | |
0a29d0f1 | 378 | // |
46d29e70 | 379 | |
380 | Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1); | |
381 | Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)- | |
382 | (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2)); | |
383 | Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)- | |
384 | (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1)); | |
385 | ||
386 | Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b); | |
387 | ||
388 | return -a/(d*y1-b)*xr/sqrt(xr*xr+yr*yr); | |
389 | } | |
390 | ||
391 | //_____________________________________________________________________ | |
5443e65e | 392 | inline Double_t f3trd(Double_t x1,Double_t y1, |
a9814c09 | 393 | Double_t x2,Double_t y2, |
394 | Double_t z1,Double_t z2) | |
46d29e70 | 395 | { |
0a29d0f1 | 396 | // |
46d29e70 | 397 | // Initial approximation of the tangent of the track dip angle |
0a29d0f1 | 398 | // |
46d29e70 | 399 | |
400 | return (z1 - z2)/sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)); | |
401 | } | |
402 | ||
46d29e70 | 403 | //___________________________________________________________________ |
5443e65e | 404 | Int_t AliTRDtracker::Clusters2Tracks(const TFile *inp, TFile *out) |
46d29e70 | 405 | { |
0a29d0f1 | 406 | // |
5443e65e | 407 | // Finds tracks within the TRD. File <inp> is expected to contain seeds |
408 | // at the outer part of the TRD. If <inp> is NULL, the seeds | |
409 | // are found within the TRD if fAddTRDseeds is TRUE. | |
410 | // The tracks are propagated to the innermost time bin | |
411 | // of the TRD and stored in file <out>. | |
0a29d0f1 | 412 | // |
a819a5f7 | 413 | |
5443e65e | 414 | LoadEvent(); |
415 | ||
416 | TDirectory *savedir=gDirectory; | |
a819a5f7 | 417 | |
5443e65e | 418 | char tname[100]; |
a819a5f7 | 419 | |
5443e65e | 420 | if (!out->IsOpen()) { |
421 | cerr<<"AliTRDtracker::Clusters2Tracks(): output file is not open !\n"; | |
422 | return 1; | |
423 | } | |
46d29e70 | 424 | |
e24ea474 | 425 | sprintf(tname,"seedTRDtoTPC_%d",GetEventNumber()); |
5443e65e | 426 | TTree tpc_tree(tname,"Tree with seeds from TRD at outer TPC pad row"); |
427 | AliTPCtrack *iotrack=0; | |
428 | tpc_tree.Branch("tracks","AliTPCtrack",&iotrack,32000,0); | |
429 | ||
e24ea474 | 430 | sprintf(tname,"TreeT%d_TRD",GetEventNumber()); |
5443e65e | 431 | TTree trd_tree(tname,"TRD tracks at inner TRD time bin"); |
432 | AliTRDtrack *iotrack_trd=0; | |
433 | trd_tree.Branch("tracks","AliTRDtrack",&iotrack_trd,32000,0); | |
434 | ||
435 | Int_t timeBins = fTrSec[0]->GetNumberOfTimeBins(); | |
436 | Float_t foundMin = fMinClustersInTrack * timeBins; | |
437 | ||
438 | if (inp) { | |
439 | TFile *in=(TFile*)inp; | |
440 | if (!in->IsOpen()) { | |
441 | cerr<<"AliTRDtracker::Clusters2Tracks(): file with seeds is not open !\n"; | |
442 | cerr<<" ... going for seeds finding inside the TRD\n"; | |
443 | } | |
444 | else { | |
445 | in->cd(); | |
e24ea474 | 446 | sprintf(tname,"TRDb_%d",GetEventNumber()); |
5443e65e | 447 | TTree *seedTree=(TTree*)in->Get(tname); |
448 | if (!seedTree) { | |
a9814c09 | 449 | cerr<<"AliTRDtracker::Clusters2Tracks(): "; |
450 | cerr<<"can't get a tree with track seeds !\n"; | |
451 | return 3; | |
5443e65e | 452 | } |
453 | AliTRDtrack *seed=new AliTRDtrack; | |
454 | seedTree->SetBranchAddress("tracks",&seed); | |
455 | ||
456 | Int_t n=(Int_t)seedTree->GetEntries(); | |
457 | for (Int_t i=0; i<n; i++) { | |
a9814c09 | 458 | seedTree->GetEvent(i); |
459 | seed->ResetCovariance(); | |
460 | AliTRDtrack *tr = new AliTRDtrack(*seed,seed->GetAlpha()); | |
461 | fSeeds->AddLast(tr); | |
462 | fNseeds++; | |
5443e65e | 463 | } |
464 | delete seed; | |
fd621f36 | 465 | delete seedTree; |
5443e65e | 466 | } |
467 | } | |
46d29e70 | 468 | |
5443e65e | 469 | out->cd(); |
46d29e70 | 470 | |
fd621f36 | 471 | |
5443e65e | 472 | // find tracks from loaded seeds |
a819a5f7 | 473 | |
5443e65e | 474 | Int_t nseed=fSeeds->GetEntriesFast(); |
475 | Int_t i, found = 0; | |
476 | Int_t innerTB = fTrSec[0]->GetInnerTimeBin(); | |
477 | ||
478 | for (i=0; i<nseed; i++) { | |
479 | AliTRDtrack *pt=(AliTRDtrack*)fSeeds->UncheckedAt(i), &t=*pt; | |
480 | FollowProlongation(t, innerTB); | |
481 | if (t.GetNumberOfClusters() >= foundMin) { | |
482 | UseClusters(&t); | |
483 | CookLabel(pt, 1-fLabelFraction); | |
484 | // t.CookdEdx(); | |
485 | } | |
486 | iotrack_trd = pt; | |
487 | trd_tree.Fill(); | |
e24ea474 | 488 | found++; |
489 | // cout<<found<<'\r'; | |
5443e65e | 490 | |
491 | if(PropagateToTPC(t)) { | |
492 | AliTPCtrack *tpc = new AliTPCtrack(*pt,pt->GetAlpha()); | |
493 | iotrack = tpc; | |
494 | tpc_tree.Fill(); | |
495 | delete tpc; | |
496 | } | |
497 | delete fSeeds->RemoveAt(i); | |
498 | fNseeds--; | |
499 | } | |
a819a5f7 | 500 | |
5443e65e | 501 | cout<<"Number of loaded seeds: "<<nseed<<endl; |
502 | cout<<"Number of found tracks from loaded seeds: "<<found<<endl; | |
a819a5f7 | 503 | |
5443e65e | 504 | // after tracks from loaded seeds are found and the corresponding |
505 | // clusters are used, look for additional seeds from TRD | |
46d29e70 | 506 | |
5443e65e | 507 | if(fAddTRDseeds) { |
508 | // Find tracks for the seeds in the TRD | |
509 | Int_t timeBins = fTrSec[0]->GetNumberOfTimeBins(); | |
510 | ||
511 | Int_t nSteps = (Int_t) (fSeedDepth / fSeedStep); | |
512 | Int_t gap = (Int_t) (timeBins * fSeedGap); | |
513 | Int_t step = (Int_t) (timeBins * fSeedStep); | |
514 | ||
515 | // make a first turn with tight cut on initial curvature | |
516 | for(Int_t turn = 1; turn <= 2; turn++) { | |
517 | if(turn == 2) { | |
a9814c09 | 518 | nSteps = (Int_t) (fSeedDepth / (3*fSeedStep)); |
519 | step = (Int_t) (timeBins * (3*fSeedStep)); | |
5443e65e | 520 | } |
521 | for(Int_t i=0; i<nSteps; i++) { | |
a9814c09 | 522 | Int_t outer=timeBins-1-i*step; |
523 | Int_t inner=outer-gap; | |
46d29e70 | 524 | |
a9814c09 | 525 | nseed=fSeeds->GetEntriesFast(); |
5443e65e | 526 | |
a9814c09 | 527 | MakeSeeds(inner, outer, turn); |
5443e65e | 528 | |
a9814c09 | 529 | nseed=fSeeds->GetEntriesFast(); |
530 | printf("\n turn %d, step %d: number of seeds for TRD inward %d\n", | |
531 | turn, i, nseed); | |
532 | ||
533 | for (Int_t i=0; i<nseed; i++) { | |
534 | AliTRDtrack *pt=(AliTRDtrack*)fSeeds->UncheckedAt(i), &t=*pt; | |
535 | FollowProlongation(t,innerTB); | |
536 | if (t.GetNumberOfClusters() >= foundMin) { | |
537 | UseClusters(&t); | |
538 | CookLabel(pt, 1-fLabelFraction); | |
539 | t.CookdEdx(); | |
e24ea474 | 540 | found++; |
541 | // cout<<found<<'\r'; | |
a9814c09 | 542 | iotrack_trd = pt; |
543 | trd_tree.Fill(); | |
544 | if(PropagateToTPC(t)) { | |
545 | AliTPCtrack *tpc = new AliTPCtrack(*pt,pt->GetAlpha()); | |
546 | iotrack = tpc; | |
547 | tpc_tree.Fill(); | |
548 | delete tpc; | |
549 | } | |
550 | } | |
551 | delete fSeeds->RemoveAt(i); | |
552 | fNseeds--; | |
553 | } | |
46d29e70 | 554 | } |
5443e65e | 555 | } |
556 | } | |
557 | tpc_tree.Write(); | |
558 | trd_tree.Write(); | |
559 | ||
560 | cout<<"Total number of found tracks: "<<found<<endl; | |
561 | ||
562 | UnloadEvent(); | |
563 | ||
564 | savedir->cd(); | |
565 | ||
566 | return 0; | |
567 | } | |
568 | ||
569 | ||
570 | ||
571 | //_____________________________________________________________________________ | |
572 | Int_t AliTRDtracker::PropagateBack(const TFile *inp, TFile *out) { | |
573 | // | |
574 | // Reads seeds from file <inp>. The seeds are AliTPCtrack's found and | |
575 | // backpropagated by the TPC tracker. Each seed is first propagated | |
576 | // to the TRD, and then its prolongation is searched in the TRD. | |
577 | // If sufficiently long continuation of the track is found in the TRD | |
578 | // the track is updated, otherwise it's stored as originaly defined | |
579 | // by the TPC tracker. | |
580 | // | |
581 | ||
5443e65e | 582 | LoadEvent(); |
583 | ||
584 | TDirectory *savedir=gDirectory; | |
585 | ||
586 | TFile *in=(TFile*)inp; | |
587 | ||
588 | if (!in->IsOpen()) { | |
589 | cerr<<"AliTRDtracker::PropagateBack(): "; | |
590 | cerr<<"file with back propagated TPC tracks is not open !\n"; | |
591 | return 1; | |
592 | } | |
593 | ||
594 | if (!out->IsOpen()) { | |
595 | cerr<<"AliTRDtracker::PropagateBack(): "; | |
596 | cerr<<"file for back propagated TRD tracks is not open !\n"; | |
597 | return 2; | |
598 | } | |
599 | ||
600 | in->cd(); | |
601 | char tname[100]; | |
e24ea474 | 602 | sprintf(tname,"seedsTPCtoTRD_%d",GetEventNumber()); |
5443e65e | 603 | TTree *seedTree=(TTree*)in->Get(tname); |
604 | if (!seedTree) { | |
605 | cerr<<"AliTRDtracker::PropagateBack(): "; | |
606 | cerr<<"can't get a tree with seeds from TPC !\n"; | |
607 | cerr<<"check if your version of TPC tracker creates tree "<<tname<<"\n"; | |
608 | return 3; | |
609 | } | |
46d29e70 | 610 | |
5443e65e | 611 | AliTPCtrack *seed=new AliTPCtrack; |
612 | seedTree->SetBranchAddress("tracks",&seed); | |
613 | ||
614 | Int_t n=(Int_t)seedTree->GetEntries(); | |
615 | for (Int_t i=0; i<n; i++) { | |
616 | seedTree->GetEvent(i); | |
617 | Int_t lbl = seed->GetLabel(); | |
618 | AliTRDtrack *tr = new AliTRDtrack(*seed,seed->GetAlpha()); | |
619 | tr->SetSeedLabel(lbl); | |
620 | fSeeds->AddLast(tr); | |
621 | fNseeds++; | |
622 | } | |
a819a5f7 | 623 | |
5443e65e | 624 | delete seed; |
625 | delete seedTree; | |
a819a5f7 | 626 | |
5443e65e | 627 | out->cd(); |
a819a5f7 | 628 | |
5443e65e | 629 | AliTPCtrack *otrack=0; |
a819a5f7 | 630 | |
e24ea474 | 631 | sprintf(tname,"seedsTRDtoTOF1_%d",GetEventNumber()); |
5443e65e | 632 | TTree tofTree1(tname,"Tracks back propagated through TPC and TRD"); |
633 | tofTree1.Branch("tracks","AliTPCtrack",&otrack,32000,0); | |
a819a5f7 | 634 | |
e24ea474 | 635 | sprintf(tname,"seedsTRDtoTOF2_%d",GetEventNumber()); |
5443e65e | 636 | TTree tofTree2(tname,"Tracks back propagated through TPC and TRD"); |
637 | tofTree2.Branch("tracks","AliTPCtrack",&otrack,32000,0); | |
46d29e70 | 638 | |
e24ea474 | 639 | sprintf(tname,"seedsTRDtoPHOS_%d",GetEventNumber()); |
5443e65e | 640 | TTree phosTree(tname,"Tracks back propagated through TPC and TRD"); |
641 | phosTree.Branch("tracks","AliTPCtrack",&otrack,32000,0); | |
a819a5f7 | 642 | |
e24ea474 | 643 | sprintf(tname,"seedsTRDtoRICH_%d",GetEventNumber()); |
5443e65e | 644 | TTree richTree(tname,"Tracks back propagated through TPC and TRD"); |
645 | richTree.Branch("tracks","AliTPCtrack",&otrack,32000,0); | |
a819a5f7 | 646 | |
e24ea474 | 647 | sprintf(tname,"TRDb_%d",GetEventNumber()); |
5443e65e | 648 | TTree trdTree(tname,"Back propagated TRD tracks at outer TRD time bin"); |
649 | AliTRDtrack *otrack_trd=0; | |
650 | trdTree.Branch("tracks","AliTRDtrack",&otrack_trd,32000,0); | |
651 | ||
9c9d2487 | 652 | if (IsStoringBarrel()) SetBarrelTree("back"); |
653 | out->cd(); | |
46d29e70 | 654 | |
9c9d2487 | 655 | Int_t found=0; |
5443e65e | 656 | Int_t nseed=fSeeds->GetEntriesFast(); |
46d29e70 | 657 | |
a9814c09 | 658 | // Float_t foundMin = fMinClustersInTrack * fTimeBinsPerPlane * fGeom->Nplan(); |
9c9d2487 | 659 | Float_t foundMin = 40; |
a819a5f7 | 660 | |
5443e65e | 661 | Int_t outermost_tb = fTrSec[0]->GetOuterTimeBin(); |
a819a5f7 | 662 | |
5443e65e | 663 | for (Int_t i=0; i<nseed; i++) { |
a819a5f7 | 664 | |
5443e65e | 665 | AliTRDtrack *ps=(AliTRDtrack*)fSeeds->UncheckedAt(i), &s=*ps; |
666 | Int_t expectedClr = FollowBackProlongation(s); | |
9c9d2487 | 667 | |
668 | if (IsStoringBarrel()) { | |
669 | StoreBarrelTrack(ps, kLastPlane, kTrackBack); | |
670 | fBarrelTree->Fill(); | |
671 | } | |
672 | ||
5443e65e | 673 | Int_t foundClr = s.GetNumberOfClusters(); |
674 | Int_t last_tb = fTrSec[0]->GetLayerNumber(s.GetX()); | |
a819a5f7 | 675 | |
fd621f36 | 676 | // printf("seed %d: found %d out of %d expected clusters, Min is %f\n", |
a9814c09 | 677 | // i, foundClr, expectedClr, foundMin); |
a819a5f7 | 678 | |
5443e65e | 679 | if (foundClr >= foundMin) { |
a9814c09 | 680 | if(foundClr >= 2) { |
681 | s.CookdEdx(); | |
682 | CookLabel(ps, 1-fLabelFraction); | |
683 | UseClusters(ps); | |
684 | } | |
0d5b5c27 | 685 | |
686 | // Propagate to outer reference plane [SR, GSI, 18.02.2003] | |
687 | ps->PropagateTo(364.8); | |
5443e65e | 688 | otrack_trd=ps; |
689 | trdTree.Fill(); | |
e24ea474 | 690 | found++; |
691 | // cout<<found<<'\r'; | |
46d29e70 | 692 | } |
a819a5f7 | 693 | |
5443e65e | 694 | if(((expectedClr < 10) && (last_tb == outermost_tb)) || |
695 | ((expectedClr >= 10) && | |
a9814c09 | 696 | (((Float_t) foundClr) / ((Float_t) expectedClr) >= |
5443e65e | 697 | fMinFractionOfFoundClusters) && (last_tb == outermost_tb))) { |
698 | ||
699 | Double_t x_tof = 375.5; | |
700 | ||
701 | if(PropagateToOuterPlane(s,x_tof)) { | |
a9814c09 | 702 | AliTPCtrack *pt = new AliTPCtrack(*ps,ps->GetAlpha()); |
703 | otrack = pt; | |
704 | tofTree1.Fill(); | |
705 | delete pt; | |
5443e65e | 706 | |
a9814c09 | 707 | x_tof = 381.5; |
5443e65e | 708 | |
a9814c09 | 709 | if(PropagateToOuterPlane(s,x_tof)) { |
710 | AliTPCtrack *pt = new AliTPCtrack(*ps,ps->GetAlpha()); | |
711 | otrack = pt; | |
712 | tofTree2.Fill(); | |
713 | delete pt; | |
714 | ||
715 | Double_t x_phos = 460.; | |
716 | ||
717 | if(PropagateToOuterPlane(s,x_phos)) { | |
718 | AliTPCtrack *pt = new AliTPCtrack(*ps,ps->GetAlpha()); | |
719 | otrack = pt; | |
720 | phosTree.Fill(); | |
721 | delete pt; | |
722 | ||
723 | Double_t x_rich = 490+1.267; | |
724 | ||
725 | if(PropagateToOuterPlane(s,x_rich)) { | |
726 | AliTPCtrack *pt = new AliTPCtrack(*ps,ps->GetAlpha()); | |
727 | otrack = pt; | |
728 | richTree.Fill(); | |
729 | delete pt; | |
730 | } | |
731 | } | |
732 | } | |
5443e65e | 733 | } |
46d29e70 | 734 | } |
735 | } | |
5443e65e | 736 | |
9c9d2487 | 737 | |
738 | out->cd(); | |
5443e65e | 739 | tofTree1.Write(); |
740 | tofTree2.Write(); | |
741 | phosTree.Write(); | |
742 | richTree.Write(); | |
743 | trdTree.Write(); | |
46d29e70 | 744 | |
9c9d2487 | 745 | if (IsStoringBarrel()) { // [SR, 03.04.2003] |
746 | fBarrelFile->cd(); | |
747 | fBarrelTree->Write(); | |
748 | fBarrelFile->Flush(); | |
749 | } | |
750 | ||
5443e65e | 751 | savedir->cd(); |
752 | cerr<<"Number of seeds: "<<nseed<<endl; | |
753 | cerr<<"Number of back propagated TRD tracks: "<<found<<endl; | |
46d29e70 | 754 | |
5443e65e | 755 | UnloadEvent(); |
46d29e70 | 756 | |
5443e65e | 757 | return 0; |
46d29e70 | 758 | |
5443e65e | 759 | } |
46d29e70 | 760 | |
79e94bf8 | 761 | //_____________________________________________________________________________ |
762 | Int_t AliTRDtracker::PropagateBack(AliESD* event) { | |
763 | // | |
764 | // Gets seeds from ESD event. The seeds are AliTPCtrack's found and | |
765 | // backpropagated by the TPC tracker. Each seed is first propagated | |
766 | // to the TRD, and then its prolongation is searched in the TRD. | |
767 | // If sufficiently long continuation of the track is found in the TRD | |
768 | // the track is updated, otherwise it's stored as originaly defined | |
769 | // by the TPC tracker. | |
770 | // | |
771 | ||
772 | Int_t found=0; | |
773 | Float_t foundMin = 40; | |
774 | ||
775 | Int_t n = event->GetNumberOfTracks(); | |
776 | for (Int_t i=0; i<n; i++) { | |
777 | AliESDtrack* seed=event->GetTrack(i); | |
778 | ULong_t status=seed->GetStatus(); | |
779 | if ( (status & AliESDtrack::kTPCout ) == 0 ) continue; | |
780 | if ( (status & AliESDtrack::kTRDout) != 0 ) continue; | |
781 | ||
782 | Int_t lbl = seed->GetLabel(); | |
783 | AliTRDtrack *track = new AliTRDtrack(*seed); | |
784 | track->SetSeedLabel(lbl); | |
785 | fNseeds++; | |
786 | ||
787 | Int_t expectedClr = FollowBackProlongation(*track); | |
788 | ||
789 | Int_t foundClr = track->GetNumberOfClusters(); | |
790 | if (foundClr >= foundMin) { | |
791 | if(foundClr >= 2) { | |
792 | track->CookdEdx(); | |
793 | // CookLabel(track, 1-fLabelFraction); | |
794 | UseClusters(track); | |
795 | } | |
796 | ||
797 | // Propagate to outer reference plane [SR, GSI, 18.02.2003] | |
798 | // track->PropagateTo(364.8); why? | |
799 | ||
800 | seed->UpdateTrackParams(track, AliESDtrack::kTRDout); | |
801 | found++; | |
802 | } | |
803 | ||
804 | } | |
805 | ||
806 | cerr<<"Number of seeds: "<<fNseeds<<endl; | |
807 | cerr<<"Number of back propagated TRD tracks: "<<found<<endl; | |
808 | ||
809 | return 0; | |
810 | ||
811 | } | |
812 | ||
bbf92647 | 813 | |
5443e65e | 814 | //--------------------------------------------------------------------------- |
815 | Int_t AliTRDtracker::FollowProlongation(AliTRDtrack& t, Int_t rf) | |
816 | { | |
817 | // Starting from current position on track=t this function tries | |
818 | // to extrapolate the track up to timeBin=0 and to confirm prolongation | |
819 | // if a close cluster is found. Returns the number of clusters | |
820 | // expected to be found in sensitive layers | |
bbf92647 | 821 | |
5443e65e | 822 | Float_t wIndex, wTB, wChi2; |
823 | Float_t wYrt, wYclosest, wYcorrect, wYwindow; | |
824 | Float_t wZrt, wZclosest, wZcorrect, wZwindow; | |
825 | Float_t wPx, wPy, wPz, wC; | |
826 | Double_t Px, Py, Pz; | |
827 | Float_t wSigmaC2, wSigmaTgl2, wSigmaY2, wSigmaZ2; | |
46d29e70 | 828 | |
5443e65e | 829 | Int_t trackIndex = t.GetLabel(); |
46d29e70 | 830 | |
5443e65e | 831 | Int_t ns=Int_t(2*TMath::Pi()/AliTRDgeometry::GetAlpha()+0.5); |
46d29e70 | 832 | |
5443e65e | 833 | Int_t try_again=fMaxGap; |
46d29e70 | 834 | |
5443e65e | 835 | Double_t alpha=t.GetAlpha(); |
9c9d2487 | 836 | TVector2::Phi_0_2pi(alpha); |
46d29e70 | 837 | |
5443e65e | 838 | Int_t s=Int_t(alpha/AliTRDgeometry::GetAlpha())%AliTRDgeometry::kNsect; |
b3a5a838 | 839 | Double_t rad_length, rho, x, dx, y, ymax, z; |
46d29e70 | 840 | |
5443e65e | 841 | Int_t expectedNumberOfClusters = 0; |
842 | Bool_t lookForCluster; | |
46d29e70 | 843 | |
5443e65e | 844 | alpha=AliTRDgeometry::GetAlpha(); // note: change in meaning |
46d29e70 | 845 | |
5443e65e | 846 | |
847 | for (Int_t nr=fTrSec[0]->GetLayerNumber(t.GetX()); nr>rf; nr--) { | |
848 | ||
849 | y = t.GetY(); z = t.GetZ(); | |
850 | ||
851 | // first propagate to the inner surface of the current time bin | |
b3a5a838 | 852 | fTrSec[s]->GetLayer(nr)->GetPropagationParameters(y,z,dx,rho,rad_length,lookForCluster); |
5443e65e | 853 | x = fTrSec[s]->GetLayer(nr)->GetX()-dx/2; y = t.GetY(); z = t.GetZ(); |
b3a5a838 | 854 | if(!t.PropagateTo(x,rad_length,rho)) break; |
5443e65e | 855 | y = t.GetY(); |
856 | ymax = x*TMath::Tan(0.5*alpha); | |
857 | if (y > ymax) { | |
858 | s = (s+1) % ns; | |
859 | if (!t.Rotate(alpha)) break; | |
b3a5a838 | 860 | if(!t.PropagateTo(x,rad_length,rho)) break; |
5443e65e | 861 | } else if (y <-ymax) { |
862 | s = (s-1+ns) % ns; | |
863 | if (!t.Rotate(-alpha)) break; | |
b3a5a838 | 864 | if(!t.PropagateTo(x,rad_length,rho)) break; |
5443e65e | 865 | } |
46d29e70 | 866 | |
5443e65e | 867 | y = t.GetY(); z = t.GetZ(); |
868 | ||
869 | // now propagate to the middle plane of the next time bin | |
b3a5a838 | 870 | fTrSec[s]->GetLayer(nr-1)->GetPropagationParameters(y,z,dx,rho,rad_length,lookForCluster); |
5443e65e | 871 | x = fTrSec[s]->GetLayer(nr-1)->GetX(); y = t.GetY(); z = t.GetZ(); |
b3a5a838 | 872 | if(!t.PropagateTo(x,rad_length,rho)) break; |
5443e65e | 873 | y = t.GetY(); |
874 | ymax = x*TMath::Tan(0.5*alpha); | |
875 | if (y > ymax) { | |
876 | s = (s+1) % ns; | |
877 | if (!t.Rotate(alpha)) break; | |
b3a5a838 | 878 | if(!t.PropagateTo(x,rad_length,rho)) break; |
5443e65e | 879 | } else if (y <-ymax) { |
880 | s = (s-1+ns) % ns; | |
881 | if (!t.Rotate(-alpha)) break; | |
b3a5a838 | 882 | if(!t.PropagateTo(x,rad_length,rho)) break; |
5443e65e | 883 | } |
46d29e70 | 884 | |
46d29e70 | 885 | |
5443e65e | 886 | if(lookForCluster) { |
a819a5f7 | 887 | |
5443e65e | 888 | expectedNumberOfClusters++; |
5443e65e | 889 | wIndex = (Float_t) t.GetLabel(); |
890 | wTB = nr; | |
46d29e70 | 891 | |
5443e65e | 892 | AliTRDpropagationLayer& time_bin=*(fTrSec[s]->GetLayer(nr-1)); |
46d29e70 | 893 | |
5443e65e | 894 | Double_t sy2=ExpectedSigmaY2(x,t.GetTgl(),t.GetPt()); |
5443e65e | 895 | Double_t sz2=ExpectedSigmaZ2(x,t.GetTgl()); |
46d29e70 | 896 | |
5443e65e | 897 | Double_t road; |
898 | if((t.GetSigmaY2() + sy2) > 0) road=10.*sqrt(t.GetSigmaY2() + sy2); | |
899 | else return expectedNumberOfClusters; | |
900 | ||
901 | wYrt = (Float_t) y; | |
902 | wZrt = (Float_t) z; | |
903 | wYwindow = (Float_t) road; | |
904 | t.GetPxPyPz(Px,Py,Pz); | |
905 | wPx = (Float_t) Px; | |
906 | wPy = (Float_t) Py; | |
907 | wPz = (Float_t) Pz; | |
908 | wC = (Float_t) t.GetC(); | |
909 | wSigmaC2 = (Float_t) t.GetSigmaC2(); | |
910 | wSigmaTgl2 = (Float_t) t.GetSigmaTgl2(); | |
911 | wSigmaY2 = (Float_t) t.GetSigmaY2(); | |
912 | wSigmaZ2 = (Float_t) t.GetSigmaZ2(); | |
913 | wChi2 = -1; | |
914 | ||
915 | if (road>fWideRoad) { | |
a9814c09 | 916 | if (t.GetNumberOfClusters()>4) |
917 | cerr<<t.GetNumberOfClusters() | |
918 | <<"FindProlongation warning: Too broad road !\n"; | |
919 | return 0; | |
5443e65e | 920 | } |
921 | ||
922 | AliTRDcluster *cl=0; | |
923 | UInt_t index=0; | |
924 | ||
925 | Double_t max_chi2=fMaxChi2; | |
926 | ||
927 | wYclosest = 12345678; | |
928 | wYcorrect = 12345678; | |
929 | wZclosest = 12345678; | |
930 | wZcorrect = 12345678; | |
931 | wZwindow = TMath::Sqrt(2.25 * 12 * sz2); | |
932 | ||
933 | // Find the closest correct cluster for debugging purposes | |
934 | if (time_bin) { | |
a9814c09 | 935 | Float_t minDY = 1000000; |
936 | for (Int_t i=0; i<time_bin; i++) { | |
937 | AliTRDcluster* c=(AliTRDcluster*)(time_bin[i]); | |
938 | if((c->GetLabel(0) != trackIndex) && | |
939 | (c->GetLabel(1) != trackIndex) && | |
940 | (c->GetLabel(2) != trackIndex)) continue; | |
941 | if(TMath::Abs(c->GetY() - y) > minDY) continue; | |
942 | minDY = TMath::Abs(c->GetY() - y); | |
943 | wYcorrect = c->GetY(); | |
944 | wZcorrect = c->GetZ(); | |
945 | ||
946 | Double_t h01 = GetTiltFactor(c); | |
947 | wChi2 = t.GetPredictedChi2(c, h01); | |
948 | } | |
5443e65e | 949 | } |
46d29e70 | 950 | |
5443e65e | 951 | // Now go for the real cluster search |
a819a5f7 | 952 | |
5443e65e | 953 | if (time_bin) { |
a819a5f7 | 954 | |
a9814c09 | 955 | for (Int_t i=time_bin.Find(y-road); i<time_bin; i++) { |
956 | AliTRDcluster* c=(AliTRDcluster*)(time_bin[i]); | |
957 | if (c->GetY() > y+road) break; | |
958 | if (c->IsUsed() > 0) continue; | |
959 | if((c->GetZ()-z)*(c->GetZ()-z) > 3 * sz2) continue; | |
960 | ||
961 | Double_t h01 = GetTiltFactor(c); | |
962 | Double_t chi2=t.GetPredictedChi2(c,h01); | |
963 | ||
964 | if (chi2 > max_chi2) continue; | |
965 | max_chi2=chi2; | |
966 | cl=c; | |
967 | index=time_bin.GetIndex(i); | |
968 | } | |
969 | ||
970 | if(!cl) { | |
971 | ||
972 | for (Int_t i=time_bin.Find(y-road); i<time_bin; i++) { | |
973 | AliTRDcluster* c=(AliTRDcluster*)(time_bin[i]); | |
974 | ||
975 | if (c->GetY() > y+road) break; | |
976 | if (c->IsUsed() > 0) continue; | |
977 | if((c->GetZ()-z)*(c->GetZ()-z) > 12 * sz2) continue; | |
978 | ||
979 | Double_t h01 = GetTiltFactor(c); | |
980 | Double_t chi2=t.GetPredictedChi2(c, h01); | |
981 | ||
982 | if (chi2 > max_chi2) continue; | |
983 | max_chi2=chi2; | |
984 | cl=c; | |
985 | index=time_bin.GetIndex(i); | |
986 | } | |
987 | } | |
988 | ||
989 | ||
990 | if (cl) { | |
991 | wYclosest = cl->GetY(); | |
992 | wZclosest = cl->GetZ(); | |
993 | Double_t h01 = GetTiltFactor(cl); | |
994 | ||
995 | t.SetSampledEdx(cl->GetQ()/dx,t.GetNumberOfClusters()); | |
996 | if(!t.Update(cl,max_chi2,index,h01)) { | |
997 | if(!try_again--) return 0; | |
998 | } | |
999 | else try_again=fMaxGap; | |
1000 | } | |
1001 | else { | |
1002 | if (try_again==0) break; | |
1003 | try_again--; | |
1004 | } | |
1005 | ||
1006 | /* | |
1007 | if((((Int_t) wTB)%15 == 0) || (((Int_t) wTB)%15 == 14)) { | |
1008 | ||
1009 | printf(" %f", wIndex); //1 | |
1010 | printf(" %f", wTB); //2 | |
1011 | printf(" %f", wYrt); //3 | |
1012 | printf(" %f", wYclosest); //4 | |
1013 | printf(" %f", wYcorrect); //5 | |
1014 | printf(" %f", wYwindow); //6 | |
1015 | printf(" %f", wZrt); //7 | |
1016 | printf(" %f", wZclosest); //8 | |
1017 | printf(" %f", wZcorrect); //9 | |
1018 | printf(" %f", wZwindow); //10 | |
1019 | printf(" %f", wPx); //11 | |
1020 | printf(" %f", wPy); //12 | |
1021 | printf(" %f", wPz); //13 | |
1022 | printf(" %f", wSigmaC2*1000000); //14 | |
1023 | printf(" %f", wSigmaTgl2*1000); //15 | |
1024 | printf(" %f", wSigmaY2); //16 | |
1025 | // printf(" %f", wSigmaZ2); //17 | |
1026 | printf(" %f", wChi2); //17 | |
1027 | printf(" %f", wC); //18 | |
1028 | printf("\n"); | |
1029 | } | |
1030 | */ | |
5443e65e | 1031 | } |
1032 | } | |
1033 | } | |
1034 | return expectedNumberOfClusters; | |
1035 | ||
1036 | ||
1037 | } | |
a819a5f7 | 1038 | |
5443e65e | 1039 | //___________________________________________________________________ |
46d29e70 | 1040 | |
5443e65e | 1041 | Int_t AliTRDtracker::FollowBackProlongation(AliTRDtrack& t) |
1042 | { | |
1043 | // Starting from current radial position of track <t> this function | |
1044 | // extrapolates the track up to outer timebin and in the sensitive | |
1045 | // layers confirms prolongation if a close cluster is found. | |
1046 | // Returns the number of clusters expected to be found in sensitive layers | |
46d29e70 | 1047 | |
5443e65e | 1048 | Float_t wIndex, wTB, wChi2; |
1049 | Float_t wYrt, wYclosest, wYcorrect, wYwindow; | |
1050 | Float_t wZrt, wZclosest, wZcorrect, wZwindow; | |
1051 | Float_t wPx, wPy, wPz, wC; | |
1052 | Double_t Px, Py, Pz; | |
1053 | Float_t wSigmaC2, wSigmaTgl2, wSigmaY2, wSigmaZ2; | |
46d29e70 | 1054 | |
5443e65e | 1055 | Int_t trackIndex = t.GetLabel(); |
5443e65e | 1056 | Int_t try_again=fMaxGap; |
46d29e70 | 1057 | |
5443e65e | 1058 | Double_t alpha=t.GetAlpha(); |
9c9d2487 | 1059 | TVector2::Phi_0_2pi(alpha); |
46d29e70 | 1060 | |
9c9d2487 | 1061 | Int_t s; |
46d29e70 | 1062 | |
5443e65e | 1063 | Int_t outerTB = fTrSec[0]->GetOuterTimeBin(); |
9c9d2487 | 1064 | Double_t rad_length, rho, x, dx, y, ymax = 0, z; |
5443e65e | 1065 | Bool_t lookForCluster; |
a819a5f7 | 1066 | |
5443e65e | 1067 | Int_t expectedNumberOfClusters = 0; |
1068 | x = t.GetX(); | |
46d29e70 | 1069 | |
5443e65e | 1070 | alpha=AliTRDgeometry::GetAlpha(); // note: change in meaning |
46d29e70 | 1071 | |
9c9d2487 | 1072 | Int_t nRefPlane = kFirstPlane; |
1073 | Bool_t isNewLayer = kFALSE; | |
46d29e70 | 1074 | |
9c9d2487 | 1075 | Double_t chi2; |
1076 | Double_t minDY; | |
46d29e70 | 1077 | |
9c9d2487 | 1078 | for (Int_t nr=fTrSec[0]->GetLayerNumber(t.GetX()); nr<outerTB+1; nr++) { |
1079 | ||
1080 | y = t.GetY(); | |
1081 | z = t.GetZ(); | |
46d29e70 | 1082 | |
5443e65e | 1083 | // first propagate to the outer surface of the current time bin |
46d29e70 | 1084 | |
9c9d2487 | 1085 | s = t.GetSector(); |
b3a5a838 | 1086 | fTrSec[s]->GetLayer(nr)->GetPropagationParameters(y,z,dx,rho,rad_length,lookForCluster); |
9c9d2487 | 1087 | x = fTrSec[s]->GetLayer(nr)->GetX()+dx/2; |
1088 | y = t.GetY(); | |
1089 | z = t.GetZ(); | |
46d29e70 | 1090 | |
b3a5a838 | 1091 | if(!t.PropagateTo(x,rad_length,rho)) break; |
9c9d2487 | 1092 | if (!AdjustSector(&t)) break; |
1093 | s = t.GetSector(); | |
1094 | if (!t.PropagateTo(x,rad_length,rho)) break; | |
46d29e70 | 1095 | |
5443e65e | 1096 | y = t.GetY(); |
9c9d2487 | 1097 | z = t.GetZ(); |
a819a5f7 | 1098 | |
9c9d2487 | 1099 | // Barrel Tracks [SR, 04.04.2003] |
46d29e70 | 1100 | |
9c9d2487 | 1101 | s = t.GetSector(); |
1102 | if (fTrSec[s]->GetLayer(nr)->IsSensitive() != | |
1103 | fTrSec[s]->GetLayer(nr+1)->IsSensitive() ) { | |
46d29e70 | 1104 | |
79e94bf8 | 1105 | // if (IsStoringBarrel()) StoreBarrelTrack(&t, nRefPlane++, kTrackBack); |
9c9d2487 | 1106 | } |
46d29e70 | 1107 | |
9c9d2487 | 1108 | if (fTrSec[s]->GetLayer(nr-1)->IsSensitive() && |
1109 | ! fTrSec[s]->GetLayer(nr)->IsSensitive()) { | |
1110 | isNewLayer = kTRUE; | |
1111 | } else {isNewLayer = kFALSE;} | |
46d29e70 | 1112 | |
5443e65e | 1113 | y = t.GetY(); |
9c9d2487 | 1114 | z = t.GetZ(); |
a819a5f7 | 1115 | |
9c9d2487 | 1116 | // now propagate to the middle plane of the next time bin |
1117 | fTrSec[s]->GetLayer(nr+1)->GetPropagationParameters(y,z,dx,rho,rad_length,lookForCluster); | |
5443e65e | 1118 | |
9c9d2487 | 1119 | x = fTrSec[s]->GetLayer(nr+1)->GetX(); |
b3a5a838 | 1120 | if(!t.PropagateTo(x,rad_length,rho)) break; |
9c9d2487 | 1121 | if (!AdjustSector(&t)) break; |
1122 | s = t.GetSector(); | |
b3a5a838 | 1123 | if(!t.PropagateTo(x,rad_length,rho)) break; |
46d29e70 | 1124 | |
9c9d2487 | 1125 | y = t.GetY(); |
1126 | z = t.GetZ(); | |
1127 | ||
1128 | if(fVocal) printf("nr+1=%d, x %f, z %f, y %f, ymax %f\n",nr+1,x,z,y,ymax); | |
5443e65e | 1129 | // printf("label %d, pl %d, lookForCluster %d \n", |
a9814c09 | 1130 | // trackIndex, nr+1, lookForCluster); |
5443e65e | 1131 | |
1132 | if(lookForCluster) { | |
1133 | expectedNumberOfClusters++; | |
1134 | ||
1135 | wIndex = (Float_t) t.GetLabel(); | |
1136 | wTB = fTrSec[s]->GetLayer(nr+1)->GetTimeBinIndex(); | |
1137 | ||
1138 | AliTRDpropagationLayer& time_bin=*(fTrSec[s]->GetLayer(nr+1)); | |
1139 | Double_t sy2=ExpectedSigmaY2(t.GetX(),t.GetTgl(),t.GetPt()); | |
1140 | Double_t sz2=ExpectedSigmaZ2(t.GetX(),t.GetTgl()); | |
1141 | if((t.GetSigmaY2() + sy2) < 0) break; | |
1142 | Double_t road = 10.*sqrt(t.GetSigmaY2() + sy2); | |
1143 | Double_t y=t.GetY(), z=t.GetZ(); | |
1144 | ||
1145 | wYrt = (Float_t) y; | |
1146 | wZrt = (Float_t) z; | |
1147 | wYwindow = (Float_t) road; | |
1148 | t.GetPxPyPz(Px,Py,Pz); | |
1149 | wPx = (Float_t) Px; | |
1150 | wPy = (Float_t) Py; | |
1151 | wPz = (Float_t) Pz; | |
1152 | wC = (Float_t) t.GetC(); | |
1153 | wSigmaC2 = (Float_t) t.GetSigmaC2(); | |
1154 | wSigmaTgl2 = (Float_t) t.GetSigmaTgl2(); | |
1155 | wSigmaY2 = (Float_t) t.GetSigmaY2(); | |
1156 | wSigmaZ2 = (Float_t) t.GetSigmaZ2(); | |
1157 | wChi2 = -1; | |
1158 | ||
1159 | if (road>fWideRoad) { | |
a9814c09 | 1160 | if (t.GetNumberOfClusters()>4) |
1161 | cerr<<t.GetNumberOfClusters() | |
1162 | <<"FindProlongation warning: Too broad road !\n"; | |
1163 | return 0; | |
5443e65e | 1164 | } |
1165 | ||
1166 | AliTRDcluster *cl=0; | |
1167 | UInt_t index=0; | |
1168 | ||
1169 | Double_t max_chi2=fMaxChi2; | |
1170 | ||
9c9d2487 | 1171 | if (isNewLayer) { |
1172 | road = 3 * road; | |
1173 | //sz2 = 3 * sz2; | |
1174 | max_chi2 = 10 * fMaxChi2; | |
1175 | } | |
1176 | ||
1177 | if (nRefPlane == kFirstPlane) max_chi2 = 20 * fMaxChi2; | |
1178 | if (nRefPlane == kFirstPlane+2) max_chi2 = 15 * fMaxChi2; | |
1179 | if (t.GetNRotate() > 0) max_chi2 = 3 * max_chi2; | |
1180 | ||
1181 | ||
5443e65e | 1182 | wYclosest = 12345678; |
1183 | wYcorrect = 12345678; | |
1184 | wZclosest = 12345678; | |
1185 | wZcorrect = 12345678; | |
1186 | wZwindow = TMath::Sqrt(2.25 * 12 * sz2); | |
1187 | ||
1188 | // Find the closest correct cluster for debugging purposes | |
1189 | if (time_bin) { | |
9c9d2487 | 1190 | minDY = 1000000; |
a9814c09 | 1191 | for (Int_t i=0; i<time_bin; i++) { |
1192 | AliTRDcluster* c=(AliTRDcluster*)(time_bin[i]); | |
1193 | if((c->GetLabel(0) != trackIndex) && | |
1194 | (c->GetLabel(1) != trackIndex) && | |
1195 | (c->GetLabel(2) != trackIndex)) continue; | |
1196 | if(TMath::Abs(c->GetY() - y) > minDY) continue; | |
9c9d2487 | 1197 | //minDY = TMath::Abs(c->GetY() - y); |
1198 | minDY = c->GetY() - y; | |
a9814c09 | 1199 | wYcorrect = c->GetY(); |
1200 | wZcorrect = c->GetZ(); | |
1201 | ||
1202 | Double_t h01 = GetTiltFactor(c); | |
1203 | wChi2 = t.GetPredictedChi2(c, h01); | |
1204 | } | |
5443e65e | 1205 | } |
46d29e70 | 1206 | |
5443e65e | 1207 | // Now go for the real cluster search |
46d29e70 | 1208 | |
5443e65e | 1209 | if (time_bin) { |
1210 | ||
a9814c09 | 1211 | for (Int_t i=time_bin.Find(y-road); i<time_bin; i++) { |
1212 | AliTRDcluster* c=(AliTRDcluster*)(time_bin[i]); | |
1213 | if (c->GetY() > y+road) break; | |
9c9d2487 | 1214 | if (c->IsUsed() > 0) continue; |
a9814c09 | 1215 | if((c->GetZ()-z)*(c->GetZ()-z) > 3 * sz2) continue; |
1216 | ||
1217 | Double_t h01 = GetTiltFactor(c); | |
9c9d2487 | 1218 | chi2=t.GetPredictedChi2(c,h01); |
a9814c09 | 1219 | |
1220 | if (chi2 > max_chi2) continue; | |
1221 | max_chi2=chi2; | |
1222 | cl=c; | |
1223 | index=time_bin.GetIndex(i); | |
9c9d2487 | 1224 | |
1225 | //check is correct | |
1226 | if((c->GetLabel(0) != trackIndex) && | |
1227 | (c->GetLabel(1) != trackIndex) && | |
1228 | (c->GetLabel(2) != trackIndex)) t.AddNWrong(); | |
a9814c09 | 1229 | } |
5443e65e | 1230 | |
a9814c09 | 1231 | if(!cl) { |
1232 | ||
1233 | for (Int_t i=time_bin.Find(y-road); i<time_bin; i++) { | |
1234 | AliTRDcluster* c=(AliTRDcluster*)(time_bin[i]); | |
1235 | ||
1236 | if (c->GetY() > y+road) break; | |
9c9d2487 | 1237 | if (c->IsUsed() > 0) continue; |
a9814c09 | 1238 | if((c->GetZ()-z)*(c->GetZ()-z) > 2.25 * 12 * sz2) continue; |
1239 | ||
1240 | Double_t h01 = GetTiltFactor(c); | |
9c9d2487 | 1241 | chi2=t.GetPredictedChi2(c,h01); |
a9814c09 | 1242 | |
1243 | if (chi2 > max_chi2) continue; | |
1244 | max_chi2=chi2; | |
1245 | cl=c; | |
1246 | index=time_bin.GetIndex(i); | |
1247 | } | |
1248 | } | |
1249 | ||
1250 | if (cl) { | |
1251 | wYclosest = cl->GetY(); | |
1252 | wZclosest = cl->GetZ(); | |
1253 | ||
1254 | t.SetSampledEdx(cl->GetQ()/dx,t.GetNumberOfClusters()); | |
1255 | Double_t h01 = GetTiltFactor(cl); | |
1256 | if(!t.Update(cl,max_chi2,index,h01)) { | |
1257 | if(!try_again--) return 0; | |
1258 | } | |
1259 | else try_again=fMaxGap; | |
1260 | } | |
1261 | else { | |
1262 | if (try_again==0) break; | |
1263 | try_again--; | |
9c9d2487 | 1264 | |
1265 | //if (minDY < 1000000 && isNewLayer) | |
1266 | //cout << "\t" << nRefPlane << "\t" << "\t" << t.GetNRotate() << "\t" << | |
1267 | // road << "\t" << minDY << "\t" << chi2 << "\t" << wChi2 << "\t" << max_chi2 << endl; | |
1268 | ||
a9814c09 | 1269 | } |
1270 | ||
9c9d2487 | 1271 | isNewLayer = kFALSE; |
1272 | ||
a9814c09 | 1273 | /* |
1274 | if((((Int_t) wTB)%15 == 0) || (((Int_t) wTB)%15 == 14)) { | |
1275 | ||
1276 | printf(" %f", wIndex); //1 | |
1277 | printf(" %f", wTB); //2 | |
1278 | printf(" %f", wYrt); //3 | |
1279 | printf(" %f", wYclosest); //4 | |
1280 | printf(" %f", wYcorrect); //5 | |
1281 | printf(" %f", wYwindow); //6 | |
1282 | printf(" %f", wZrt); //7 | |
1283 | printf(" %f", wZclosest); //8 | |
1284 | printf(" %f", wZcorrect); //9 | |
1285 | printf(" %f", wZwindow); //10 | |
1286 | printf(" %f", wPx); //11 | |
1287 | printf(" %f", wPy); //12 | |
1288 | printf(" %f", wPz); //13 | |
1289 | printf(" %f", wSigmaC2*1000000); //14 | |
1290 | printf(" %f", wSigmaTgl2*1000); //15 | |
1291 | printf(" %f", wSigmaY2); //16 | |
1292 | // printf(" %f", wSigmaZ2); //17 | |
1293 | printf(" %f", wChi2); //17 | |
1294 | printf(" %f", wC); //18 | |
1295 | printf("\n"); | |
1296 | } | |
1297 | */ | |
5443e65e | 1298 | } |
1299 | } | |
1300 | } | |
1301 | return expectedNumberOfClusters; | |
9c9d2487 | 1302 | |
1303 | ||
5443e65e | 1304 | } |
1305 | ||
1306 | //___________________________________________________________________ | |
1307 | ||
1308 | Int_t AliTRDtracker::PropagateToOuterPlane(AliTRDtrack& t, Double_t xToGo) | |
1309 | { | |
1310 | // Starting from current radial position of track <t> this function | |
1311 | // extrapolates the track up to radial position <xToGo>. | |
1312 | // Returns 1 if track reaches the plane, and 0 otherwise | |
1313 | ||
1314 | Int_t ns=Int_t(2*TMath::Pi()/AliTRDgeometry::GetAlpha()+0.5); | |
1315 | ||
1316 | Double_t alpha=t.GetAlpha(); | |
1317 | ||
1318 | if (alpha > 2.*TMath::Pi()) alpha -= 2.*TMath::Pi(); | |
1319 | if (alpha < 0. ) alpha += 2.*TMath::Pi(); | |
1320 | ||
1321 | Int_t s=Int_t(alpha/AliTRDgeometry::GetAlpha())%AliTRDgeometry::kNsect; | |
1322 | ||
1323 | Bool_t lookForCluster; | |
b3a5a838 | 1324 | Double_t rad_length, rho, x, dx, y, ymax, z; |
5443e65e | 1325 | |
1326 | x = t.GetX(); | |
1327 | ||
1328 | alpha=AliTRDgeometry::GetAlpha(); // note: change in meaning | |
1329 | ||
1330 | Int_t plToGo = fTrSec[0]->GetLayerNumber(xToGo); | |
1331 | ||
1332 | for (Int_t nr=fTrSec[0]->GetLayerNumber(x); nr<plToGo; nr++) { | |
1333 | ||
1334 | y = t.GetY(); z = t.GetZ(); | |
1335 | ||
1336 | // first propagate to the outer surface of the current time bin | |
b3a5a838 | 1337 | fTrSec[s]->GetLayer(nr)->GetPropagationParameters(y,z,dx,rho,rad_length,lookForCluster); |
5443e65e | 1338 | x = fTrSec[s]->GetLayer(nr)->GetX()+dx/2; y = t.GetY(); z = t.GetZ(); |
b3a5a838 | 1339 | if(!t.PropagateTo(x,rad_length,rho)) return 0; |
5443e65e | 1340 | y = t.GetY(); |
1341 | ymax = x*TMath::Tan(0.5*alpha); | |
1342 | if (y > ymax) { | |
1343 | s = (s+1) % ns; | |
1344 | if (!t.Rotate(alpha)) return 0; | |
1345 | } else if (y <-ymax) { | |
1346 | s = (s-1+ns) % ns; | |
1347 | if (!t.Rotate(-alpha)) return 0; | |
1348 | } | |
b3a5a838 | 1349 | if(!t.PropagateTo(x,rad_length,rho)) return 0; |
5443e65e | 1350 | |
1351 | y = t.GetY(); z = t.GetZ(); | |
1352 | ||
1353 | // now propagate to the middle plane of the next time bin | |
b3a5a838 | 1354 | fTrSec[s]->GetLayer(nr+1)->GetPropagationParameters(y,z,dx,rho,rad_length,lookForCluster); |
5443e65e | 1355 | x = fTrSec[s]->GetLayer(nr+1)->GetX(); y = t.GetY(); z = t.GetZ(); |
b3a5a838 | 1356 | if(!t.PropagateTo(x,rad_length,rho)) return 0; |
5443e65e | 1357 | y = t.GetY(); |
1358 | ymax = x*TMath::Tan(0.5*alpha); | |
1359 | if (y > ymax) { | |
1360 | s = (s+1) % ns; | |
1361 | if (!t.Rotate(alpha)) return 0; | |
1362 | } else if (y <-ymax) { | |
1363 | s = (s-1+ns) % ns; | |
1364 | if (!t.Rotate(-alpha)) return 0; | |
1365 | } | |
b3a5a838 | 1366 | if(!t.PropagateTo(x,rad_length,rho)) return 0; |
5443e65e | 1367 | } |
1368 | return 1; | |
1369 | } | |
1370 | ||
1371 | //___________________________________________________________________ | |
1372 | ||
1373 | Int_t AliTRDtracker::PropagateToTPC(AliTRDtrack& t) | |
1374 | { | |
1375 | // Starting from current radial position of track <t> this function | |
1376 | // extrapolates the track up to radial position of the outermost | |
1377 | // padrow of the TPC. | |
1378 | // Returns 1 if track reaches the TPC, and 0 otherwise | |
1379 | ||
1380 | Int_t ns=Int_t(2*TMath::Pi()/AliTRDgeometry::GetAlpha()+0.5); | |
1381 | ||
1382 | Double_t alpha=t.GetAlpha(); | |
9c9d2487 | 1383 | TVector2::Phi_0_2pi(alpha); |
5443e65e | 1384 | |
1385 | Int_t s=Int_t(alpha/AliTRDgeometry::GetAlpha())%AliTRDgeometry::kNsect; | |
1386 | ||
1387 | Bool_t lookForCluster; | |
b3a5a838 | 1388 | Double_t rad_length, rho, x, dx, y, ymax, z; |
5443e65e | 1389 | |
1390 | x = t.GetX(); | |
1391 | ||
1392 | alpha=AliTRDgeometry::GetAlpha(); // note: change in meaning | |
5443e65e | 1393 | Int_t plTPC = fTrSec[0]->GetLayerNumber(246.055); |
1394 | ||
1395 | for (Int_t nr=fTrSec[0]->GetLayerNumber(x); nr>plTPC; nr--) { | |
1396 | ||
9c9d2487 | 1397 | y = t.GetY(); |
1398 | z = t.GetZ(); | |
5443e65e | 1399 | |
1400 | // first propagate to the outer surface of the current time bin | |
b3a5a838 | 1401 | fTrSec[s]->GetLayer(nr)->GetPropagationParameters(y,z,dx,rho,rad_length,lookForCluster); |
9c9d2487 | 1402 | x = fTrSec[s]->GetLayer(nr)->GetX()-dx/2; |
1403 | ||
b3a5a838 | 1404 | if(!t.PropagateTo(x,rad_length,rho)) return 0; |
9c9d2487 | 1405 | AdjustSector(&t); |
b3a5a838 | 1406 | if(!t.PropagateTo(x,rad_length,rho)) return 0; |
5443e65e | 1407 | |
9c9d2487 | 1408 | y = t.GetY(); |
1409 | z = t.GetZ(); | |
5443e65e | 1410 | |
1411 | // now propagate to the middle plane of the next time bin | |
b3a5a838 | 1412 | fTrSec[s]->GetLayer(nr-1)->GetPropagationParameters(y,z,dx,rho,rad_length,lookForCluster); |
9c9d2487 | 1413 | x = fTrSec[s]->GetLayer(nr-1)->GetX(); |
1414 | ||
b3a5a838 | 1415 | if(!t.PropagateTo(x,rad_length,rho)) return 0; |
9c9d2487 | 1416 | AdjustSector(&t); |
b3a5a838 | 1417 | if(!t.PropagateTo(x,rad_length,rho)) return 0; |
5443e65e | 1418 | } |
1419 | return 1; | |
1420 | } | |
1421 | ||
5443e65e | 1422 | //_____________________________________________________________________________ |
1423 | void AliTRDtracker::LoadEvent() | |
1424 | { | |
1425 | // Fills clusters into TRD tracking_sectors | |
1426 | // Note that the numbering scheme for the TRD tracking_sectors | |
1427 | // differs from that of TRD sectors | |
1428 | ||
1429 | ReadClusters(fClusters); | |
1430 | Int_t ncl=fClusters->GetEntriesFast(); | |
1431 | cout<<"\n LoadSectors: sorting "<<ncl<<" clusters"<<endl; | |
1432 | ||
1433 | UInt_t index; | |
1434 | while (ncl--) { | |
e24ea474 | 1435 | // printf("\r %d left ",ncl); |
5443e65e | 1436 | AliTRDcluster *c=(AliTRDcluster*)fClusters->UncheckedAt(ncl); |
1437 | Int_t detector=c->GetDetector(), local_time_bin=c->GetLocalTimeBin(); | |
1438 | Int_t sector=fGeom->GetSector(detector); | |
1439 | Int_t plane=fGeom->GetPlane(detector); | |
1440 | ||
1441 | Int_t tracking_sector = CookSectorIndex(sector); | |
1442 | ||
1443 | Int_t gtb = fTrSec[tracking_sector]->CookTimeBinIndex(plane,local_time_bin); | |
b3a5a838 | 1444 | if(gtb < 0) continue; |
5443e65e | 1445 | Int_t layer = fTrSec[tracking_sector]->GetLayerNumber(gtb); |
1446 | ||
1447 | index=ncl; | |
1448 | fTrSec[tracking_sector]->GetLayer(layer)->InsertCluster(c,index); | |
1449 | } | |
1450 | printf("\r\n"); | |
1451 | ||
1452 | } | |
1453 | ||
1454 | //_____________________________________________________________________________ | |
1455 | void AliTRDtracker::UnloadEvent() | |
1456 | { | |
1457 | // | |
1458 | // Clears the arrays of clusters and tracks. Resets sectors and timebins | |
1459 | // | |
1460 | ||
1461 | Int_t i, nentr; | |
1462 | ||
1463 | nentr = fClusters->GetEntriesFast(); | |
1464 | for (i = 0; i < nentr; i++) delete fClusters->RemoveAt(i); | |
1465 | ||
1466 | nentr = fSeeds->GetEntriesFast(); | |
1467 | for (i = 0; i < nentr; i++) delete fSeeds->RemoveAt(i); | |
1468 | ||
1469 | nentr = fTracks->GetEntriesFast(); | |
1470 | for (i = 0; i < nentr; i++) delete fTracks->RemoveAt(i); | |
1471 | ||
1472 | Int_t nsec = AliTRDgeometry::kNsect; | |
1473 | ||
1474 | for (i = 0; i < nsec; i++) { | |
1475 | for(Int_t pl = 0; pl < fTrSec[i]->GetNumberOfLayers(); pl++) { | |
1476 | fTrSec[i]->GetLayer(pl)->Clear(); | |
1477 | } | |
1478 | } | |
1479 | ||
1480 | } | |
1481 | ||
1482 | //__________________________________________________________________________ | |
1483 | void AliTRDtracker::MakeSeeds(Int_t inner, Int_t outer, Int_t turn) | |
1484 | { | |
1485 | // Creates track seeds using clusters in timeBins=i1,i2 | |
1486 | ||
1487 | if(turn > 2) { | |
1488 | cerr<<"MakeSeeds: turn "<<turn<<" exceeds the limit of 2"<<endl; | |
1489 | return; | |
1490 | } | |
1491 | ||
1492 | Double_t x[5], c[15]; | |
1493 | Int_t max_sec=AliTRDgeometry::kNsect; | |
1494 | ||
1495 | Double_t alpha=AliTRDgeometry::GetAlpha(); | |
1496 | Double_t shift=AliTRDgeometry::GetAlpha()/2.; | |
1497 | Double_t cs=cos(alpha), sn=sin(alpha); | |
1498 | Double_t cs2=cos(2.*alpha), sn2=sin(2.*alpha); | |
1499 | ||
1500 | ||
1501 | Int_t i2 = fTrSec[0]->GetLayerNumber(inner); | |
1502 | Int_t i1 = fTrSec[0]->GetLayerNumber(outer); | |
1503 | ||
1504 | Double_t x1 =fTrSec[0]->GetX(i1); | |
1505 | Double_t xx2=fTrSec[0]->GetX(i2); | |
1506 | ||
1507 | for (Int_t ns=0; ns<max_sec; ns++) { | |
1508 | ||
1509 | Int_t nl2 = *(fTrSec[(ns-2+max_sec)%max_sec]->GetLayer(i2)); | |
1510 | Int_t nl=(*fTrSec[(ns-1+max_sec)%max_sec]->GetLayer(i2)); | |
1511 | Int_t nm=(*fTrSec[ns]->GetLayer(i2)); | |
1512 | Int_t nu=(*fTrSec[(ns+1)%max_sec]->GetLayer(i2)); | |
1513 | Int_t nu2=(*fTrSec[(ns+2)%max_sec]->GetLayer(i2)); | |
1514 | ||
1515 | AliTRDpropagationLayer& r1=*(fTrSec[ns]->GetLayer(i1)); | |
1516 | ||
1517 | for (Int_t is=0; is < r1; is++) { | |
1518 | Double_t y1=r1[is]->GetY(), z1=r1[is]->GetZ(); | |
1519 | ||
1520 | for (Int_t js=0; js < nl2+nl+nm+nu+nu2; js++) { | |
a9814c09 | 1521 | |
1522 | const AliTRDcluster *cl; | |
1523 | Double_t x2, y2, z2; | |
1524 | Double_t x3=0., y3=0.; | |
1525 | ||
1526 | if (js<nl2) { | |
1527 | if(turn != 2) continue; | |
1528 | AliTRDpropagationLayer& r2=*(fTrSec[(ns-2+max_sec)%max_sec]->GetLayer(i2)); | |
1529 | cl=r2[js]; | |
1530 | y2=cl->GetY(); z2=cl->GetZ(); | |
1531 | ||
1532 | x2= xx2*cs2+y2*sn2; | |
1533 | y2=-xx2*sn2+y2*cs2; | |
1534 | } | |
1535 | else if (js<nl2+nl) { | |
1536 | if(turn != 1) continue; | |
1537 | AliTRDpropagationLayer& r2=*(fTrSec[(ns-1+max_sec)%max_sec]->GetLayer(i2)); | |
1538 | cl=r2[js-nl2]; | |
1539 | y2=cl->GetY(); z2=cl->GetZ(); | |
1540 | ||
1541 | x2= xx2*cs+y2*sn; | |
1542 | y2=-xx2*sn+y2*cs; | |
1543 | } | |
1544 | else if (js<nl2+nl+nm) { | |
1545 | if(turn != 1) continue; | |
1546 | AliTRDpropagationLayer& r2=*(fTrSec[ns]->GetLayer(i2)); | |
1547 | cl=r2[js-nl2-nl]; | |
1548 | x2=xx2; y2=cl->GetY(); z2=cl->GetZ(); | |
1549 | } | |
1550 | else if (js<nl2+nl+nm+nu) { | |
1551 | if(turn != 1) continue; | |
1552 | AliTRDpropagationLayer& r2=*(fTrSec[(ns+1)%max_sec]->GetLayer(i2)); | |
1553 | cl=r2[js-nl2-nl-nm]; | |
1554 | y2=cl->GetY(); z2=cl->GetZ(); | |
1555 | ||
1556 | x2=xx2*cs-y2*sn; | |
1557 | y2=xx2*sn+y2*cs; | |
1558 | } | |
1559 | else { | |
1560 | if(turn != 2) continue; | |
1561 | AliTRDpropagationLayer& r2=*(fTrSec[(ns+2)%max_sec]->GetLayer(i2)); | |
1562 | cl=r2[js-nl2-nl-nm-nu]; | |
1563 | y2=cl->GetY(); z2=cl->GetZ(); | |
1564 | ||
1565 | x2=xx2*cs2-y2*sn2; | |
1566 | y2=xx2*sn2+y2*cs2; | |
1567 | } | |
1568 | ||
1569 | if(TMath::Abs(z1-z2) > fMaxSeedDeltaZ12) continue; | |
1570 | ||
1571 | Double_t zz=z1 - z1/x1*(x1-x2); | |
1572 | ||
1573 | if (TMath::Abs(zz-z2)>fMaxSeedDeltaZ) continue; | |
1574 | ||
1575 | Double_t d=(x2-x1)*(0.-y2)-(0.-x2)*(y2-y1); | |
1576 | if (d==0.) {cerr<<"TRD MakeSeeds: Straight seed !\n"; continue;} | |
1577 | ||
1578 | x[0]=y1; | |
1579 | x[1]=z1; | |
1580 | x[4]=f1trd(x1,y1,x2,y2,x3,y3); | |
1581 | ||
1582 | if (TMath::Abs(x[4]) > fMaxSeedC) continue; | |
1583 | ||
1584 | x[2]=f2trd(x1,y1,x2,y2,x3,y3); | |
1585 | ||
1586 | if (TMath::Abs(x[4]*x1-x[2]) >= 0.99999) continue; | |
1587 | ||
1588 | x[3]=f3trd(x1,y1,x2,y2,z1,z2); | |
1589 | ||
1590 | if (TMath::Abs(x[3]) > fMaxSeedTan) continue; | |
1591 | ||
1592 | Double_t a=asin(x[2]); | |
1593 | Double_t zv=z1 - x[3]/x[4]*(a+asin(x[4]*x1-x[2])); | |
1594 | ||
1595 | if (TMath::Abs(zv)>fMaxSeedVertexZ) continue; | |
1596 | ||
1597 | Double_t sy1=r1[is]->GetSigmaY2(), sz1=r1[is]->GetSigmaZ2(); | |
1598 | Double_t sy2=cl->GetSigmaY2(), sz2=cl->GetSigmaZ2(); | |
1599 | Double_t sy3=fSeedErrorSY3, sy=fSeedErrorSY, sz=fSeedErrorSZ; | |
1600 | ||
1601 | // Tilt changes | |
1602 | Double_t h01 = GetTiltFactor(r1[is]); | |
b8dc2353 | 1603 | Double_t xu_factor = 100.; |
1604 | if(fNoTilt) { | |
1605 | h01 = 0; | |
1606 | xu_factor = 1; | |
1607 | } | |
1608 | ||
fd621f36 | 1609 | sy1=sy1+sz1*h01*h01; |
b3a5a838 | 1610 | Double_t syz=sz1*(-h01); |
a9814c09 | 1611 | // end of tilt changes |
1612 | ||
b3a5a838 | 1613 | Double_t f40=(f1trd(x1,y1+sy,x2,y2,x3,y3)-x[4])/sy; |
1614 | Double_t f42=(f1trd(x1,y1,x2,y2+sy,x3,y3)-x[4])/sy; | |
1615 | Double_t f43=(f1trd(x1,y1,x2,y2,x3,y3+sy)-x[4])/sy; | |
1616 | Double_t f20=(f2trd(x1,y1+sy,x2,y2,x3,y3)-x[2])/sy; | |
1617 | Double_t f22=(f2trd(x1,y1,x2,y2+sy,x3,y3)-x[2])/sy; | |
1618 | Double_t f23=(f2trd(x1,y1,x2,y2,x3,y3+sy)-x[2])/sy; | |
1619 | Double_t f30=(f3trd(x1,y1+sy,x2,y2,z1,z2)-x[3])/sy; | |
1620 | Double_t f31=(f3trd(x1,y1,x2,y2,z1+sz,z2)-x[3])/sz; | |
1621 | Double_t f32=(f3trd(x1,y1,x2,y2+sy,z1,z2)-x[3])/sy; | |
1622 | Double_t f34=(f3trd(x1,y1,x2,y2,z1,z2+sz)-x[3])/sz; | |
1623 | ||
a9814c09 | 1624 | |
b3a5a838 | 1625 | c[0]=sy1; |
a9814c09 | 1626 | // c[1]=0.; c[2]=sz1; |
b8dc2353 | 1627 | c[1]=syz; c[2]=sz1*xu_factor; |
b3a5a838 | 1628 | c[3]=f20*sy1; c[4]=0.; c[5]=f20*sy1*f20+f22*sy2*f22+f23*sy3*f23; |
1629 | c[6]=f30*sy1; c[7]=f31*sz1; c[8]=f30*sy1*f20+f32*sy2*f22; | |
1630 | c[9]=f30*sy1*f30+f31*sz1*f31+f32*sy2*f32+f34*sz2*f34; | |
1631 | c[10]=f40*sy1; c[11]=0.; c[12]=f40*sy1*f20+f42*sy2*f22+f43*sy3*f23; | |
1632 | c[13]=f30*sy1*f40+f32*sy2*f42; | |
1633 | c[14]=f40*sy1*f40+f42*sy2*f42+f43*sy3*f43; | |
a9814c09 | 1634 | |
1635 | UInt_t index=r1.GetIndex(is); | |
1636 | ||
1637 | AliTRDtrack *track=new AliTRDtrack(r1[is],index,x,c,x1,ns*alpha+shift); | |
1638 | ||
1639 | Int_t rc=FollowProlongation(*track, i2); | |
1640 | ||
1641 | if ((rc < 1) || | |
1642 | (track->GetNumberOfClusters() < | |
1643 | (outer-inner)*fMinClustersInSeed)) delete track; | |
1644 | else { | |
1645 | fSeeds->AddLast(track); fNseeds++; | |
e24ea474 | 1646 | // cerr<<"\r found seed "<<fNseeds; |
a9814c09 | 1647 | } |
5443e65e | 1648 | } |
1649 | } | |
1650 | } | |
1651 | } | |
1652 | ||
1653 | //_____________________________________________________________________________ | |
1654 | void AliTRDtracker::ReadClusters(TObjArray *array, const TFile *inp) | |
1655 | { | |
1656 | // | |
a819a5f7 | 1657 | // Reads AliTRDclusters (option >= 0) or AliTRDrecPoints (option < 0) |
1658 | // from the file. The names of the cluster tree and branches | |
1659 | // should match the ones used in AliTRDclusterizer::WriteClusters() | |
1660 | // | |
46d29e70 | 1661 | |
a819a5f7 | 1662 | TDirectory *savedir=gDirectory; |
1663 | ||
5443e65e | 1664 | if (inp) { |
1665 | TFile *in=(TFile*)inp; | |
1666 | if (!in->IsOpen()) { | |
1667 | cerr<<"AliTRDtracker::ReadClusters(): input file is not open !\n"; | |
1668 | return; | |
1669 | } | |
1670 | else{ | |
1671 | in->cd(); | |
1672 | } | |
1673 | } | |
a819a5f7 | 1674 | |
abaf1f1d | 1675 | Char_t treeName[12]; |
e24ea474 | 1676 | sprintf(treeName,"TreeR%d_TRD",GetEventNumber()); |
5443e65e | 1677 | TTree *ClusterTree = (TTree*) gDirectory->Get(treeName); |
1678 | ||
1679 | TObjArray *ClusterArray = new TObjArray(400); | |
1680 | ||
1681 | ClusterTree->GetBranch("TRDcluster")->SetAddress(&ClusterArray); | |
1682 | ||
1683 | Int_t nEntries = (Int_t) ClusterTree->GetEntries(); | |
1684 | printf("found %d entries in %s.\n",nEntries,ClusterTree->GetName()); | |
a819a5f7 | 1685 | |
a819a5f7 | 1686 | // Loop through all entries in the tree |
1687 | Int_t nbytes; | |
1688 | AliTRDcluster *c = 0; | |
5443e65e | 1689 | printf("\n"); |
a819a5f7 | 1690 | |
1691 | for (Int_t iEntry = 0; iEntry < nEntries; iEntry++) { | |
1692 | ||
1693 | // Import the tree | |
5443e65e | 1694 | nbytes += ClusterTree->GetEvent(iEntry); |
1695 | ||
a819a5f7 | 1696 | // Get the number of points in the detector |
5443e65e | 1697 | Int_t nCluster = ClusterArray->GetEntriesFast(); |
e24ea474 | 1698 | // printf("\r Read %d clusters from entry %d", nCluster, iEntry); |
5443e65e | 1699 | |
a819a5f7 | 1700 | // Loop through all TRD digits |
1701 | for (Int_t iCluster = 0; iCluster < nCluster; iCluster++) { | |
5443e65e | 1702 | c = (AliTRDcluster*)ClusterArray->UncheckedAt(iCluster); |
a819a5f7 | 1703 | AliTRDcluster *co = new AliTRDcluster(*c); |
1704 | co->SetSigmaY2(c->GetSigmaY2() * fSY2corr); | |
5443e65e | 1705 | Int_t ltb = co->GetLocalTimeBin(); |
b8dc2353 | 1706 | if(ltb == 19) co->SetSigmaZ2(c->GetSigmaZ2()); |
1707 | else if(fNoTilt) co->SetSigmaZ2(c->GetSigmaZ2() * fSZ2corr); | |
a819a5f7 | 1708 | array->AddLast(co); |
5443e65e | 1709 | delete ClusterArray->RemoveAt(iCluster); |
a819a5f7 | 1710 | } |
1711 | } | |
1712 | ||
5443e65e | 1713 | delete ClusterArray; |
1714 | savedir->cd(); | |
1715 | ||
a819a5f7 | 1716 | } |
1717 | ||
5443e65e | 1718 | //______________________________________________________________________ |
1719 | void AliTRDtracker::ReadClusters(TObjArray *array, const Char_t *filename) | |
a819a5f7 | 1720 | { |
1721 | // | |
5443e65e | 1722 | // Reads AliTRDclusters from file <filename>. The names of the cluster |
1723 | // tree and branches should match the ones used in | |
1724 | // AliTRDclusterizer::WriteClusters() | |
1725 | // if <array> == 0, clusters are added into AliTRDtracker fCluster array | |
a819a5f7 | 1726 | // |
1727 | ||
5443e65e | 1728 | TDirectory *savedir=gDirectory; |
46d29e70 | 1729 | |
5443e65e | 1730 | TFile *file = TFile::Open(filename); |
1731 | if (!file->IsOpen()) { | |
1732 | cerr<<"Can't open file with TRD clusters"<<endl; | |
1733 | return; | |
1734 | } | |
46d29e70 | 1735 | |
5443e65e | 1736 | Char_t treeName[12]; |
e24ea474 | 1737 | sprintf(treeName,"TreeR%d_TRD",GetEventNumber()); |
5443e65e | 1738 | TTree *ClusterTree = (TTree*) gDirectory->Get(treeName); |
46d29e70 | 1739 | |
5443e65e | 1740 | if (!ClusterTree) { |
1741 | cerr<<"AliTRDtracker::ReadClusters(): "; | |
1742 | cerr<<"can't get a tree with clusters !\n"; | |
1743 | return; | |
1744 | } | |
46d29e70 | 1745 | |
5443e65e | 1746 | TObjArray *ClusterArray = new TObjArray(400); |
46d29e70 | 1747 | |
5443e65e | 1748 | ClusterTree->GetBranch("TRDcluster")->SetAddress(&ClusterArray); |
46d29e70 | 1749 | |
5443e65e | 1750 | Int_t nEntries = (Int_t) ClusterTree->GetEntries(); |
1751 | cout<<"found "<<nEntries<<" in ClusterTree"<<endl; | |
a819a5f7 | 1752 | |
5443e65e | 1753 | // Loop through all entries in the tree |
1754 | Int_t nbytes; | |
1755 | AliTRDcluster *c = 0; | |
bbf92647 | 1756 | |
5443e65e | 1757 | printf("\n"); |
bbf92647 | 1758 | |
5443e65e | 1759 | for (Int_t iEntry = 0; iEntry < nEntries; iEntry++) { |
46d29e70 | 1760 | |
5443e65e | 1761 | // Import the tree |
1762 | nbytes += ClusterTree->GetEvent(iEntry); | |
0a29d0f1 | 1763 | |
5443e65e | 1764 | // Get the number of points in the detector |
1765 | Int_t nCluster = ClusterArray->GetEntriesFast(); | |
b3a5a838 | 1766 | printf("\n Read %d clusters from entry %d", nCluster, iEntry); |
5443e65e | 1767 | |
1768 | // Loop through all TRD digits | |
1769 | for (Int_t iCluster = 0; iCluster < nCluster; iCluster++) { | |
1770 | c = (AliTRDcluster*)ClusterArray->UncheckedAt(iCluster); | |
1771 | AliTRDcluster *co = new AliTRDcluster(*c); | |
1772 | co->SetSigmaY2(c->GetSigmaY2() * fSY2corr); | |
1773 | Int_t ltb = co->GetLocalTimeBin(); | |
b8dc2353 | 1774 | if(ltb == 19) co->SetSigmaZ2(c->GetSigmaZ2()); |
1775 | else if(fNoTilt) co->SetSigmaZ2(c->GetSigmaZ2() * fSZ2corr); | |
5443e65e | 1776 | array->AddLast(co); |
1777 | delete ClusterArray->RemoveAt(iCluster); | |
1778 | } | |
46d29e70 | 1779 | } |
0a29d0f1 | 1780 | |
5443e65e | 1781 | file->Close(); |
1782 | delete ClusterArray; | |
1783 | savedir->cd(); | |
1784 | ||
1785 | } | |
1786 | ||
46d29e70 | 1787 | |
1788 | //__________________________________________________________________ | |
5443e65e | 1789 | void AliTRDtracker::CookLabel(AliKalmanTrack* pt, Float_t wrong) const { |
46d29e70 | 1790 | |
1791 | Int_t label=123456789, index, i, j; | |
5443e65e | 1792 | Int_t ncl=pt->GetNumberOfClusters(); |
1793 | const Int_t range = fTrSec[0]->GetOuterTimeBin()+1; | |
1794 | ||
1795 | Bool_t label_added; | |
46d29e70 | 1796 | |
5443e65e | 1797 | // Int_t s[range][2]; |
1798 | Int_t **s = new Int_t* [range]; | |
1799 | for (i=0; i<range; i++) { | |
d1b06c24 | 1800 | s[i] = new Int_t[2]; |
1801 | } | |
5443e65e | 1802 | for (i=0; i<range; i++) { |
46d29e70 | 1803 | s[i][0]=-1; |
1804 | s[i][1]=0; | |
1805 | } | |
1806 | ||
1807 | Int_t t0,t1,t2; | |
1808 | for (i=0; i<ncl; i++) { | |
5443e65e | 1809 | index=pt->GetClusterIndex(i); |
46d29e70 | 1810 | AliTRDcluster *c=(AliTRDcluster*)fClusters->UncheckedAt(index); |
5443e65e | 1811 | t0=c->GetLabel(0); |
1812 | t1=c->GetLabel(1); | |
1813 | t2=c->GetLabel(2); | |
46d29e70 | 1814 | } |
1815 | ||
1816 | for (i=0; i<ncl; i++) { | |
5443e65e | 1817 | index=pt->GetClusterIndex(i); |
46d29e70 | 1818 | AliTRDcluster *c=(AliTRDcluster*)fClusters->UncheckedAt(index); |
1819 | for (Int_t k=0; k<3; k++) { | |
5443e65e | 1820 | label=c->GetLabel(k); |
1821 | label_added=kFALSE; j=0; | |
46d29e70 | 1822 | if (label >= 0) { |
a9814c09 | 1823 | while ( (!label_added) && ( j < range ) ) { |
1824 | if (s[j][0]==label || s[j][1]==0) { | |
1825 | s[j][0]=label; | |
1826 | s[j][1]=s[j][1]+1; | |
1827 | label_added=kTRUE; | |
1828 | } | |
1829 | j++; | |
1830 | } | |
46d29e70 | 1831 | } |
1832 | } | |
1833 | } | |
1834 | ||
46d29e70 | 1835 | Int_t max=0; |
1836 | label = -123456789; | |
1837 | ||
5443e65e | 1838 | for (i=0; i<range; i++) { |
46d29e70 | 1839 | if (s[i][1]>max) { |
1840 | max=s[i][1]; label=s[i][0]; | |
1841 | } | |
1842 | } | |
5443e65e | 1843 | |
1844 | for (i=0; i<range; i++) { | |
1845 | delete []s[i]; | |
1846 | } | |
1847 | ||
d1b06c24 | 1848 | delete []s; |
5443e65e | 1849 | |
1850 | if ((1.- Float_t(max)/ncl) > wrong) label=-label; | |
1851 | ||
1852 | pt->SetLabel(label); | |
1853 | ||
46d29e70 | 1854 | } |
1855 | ||
5443e65e | 1856 | |
1857 | //__________________________________________________________________ | |
ca6c93d6 | 1858 | void AliTRDtracker::UseClusters(const AliKalmanTrack* t, Int_t from) const { |
5443e65e | 1859 | Int_t ncl=t->GetNumberOfClusters(); |
1860 | for (Int_t i=from; i<ncl; i++) { | |
1861 | Int_t index = t->GetClusterIndex(i); | |
1862 | AliTRDcluster *c=(AliTRDcluster*)fClusters->UncheckedAt(index); | |
1863 | c->Use(); | |
1864 | } | |
1865 | } | |
1866 | ||
1867 | ||
1868 | //_____________________________________________________________________ | |
1869 | Double_t AliTRDtracker::ExpectedSigmaY2(Double_t r, Double_t tgl, Double_t pt) | |
1870 | { | |
1871 | // Parametrised "expected" error of the cluster reconstruction in Y | |
1872 | ||
1873 | Double_t s = 0.08 * 0.08; | |
1874 | return s; | |
1875 | } | |
1876 | ||
1877 | //_____________________________________________________________________ | |
1878 | Double_t AliTRDtracker::ExpectedSigmaZ2(Double_t r, Double_t tgl) | |
0a29d0f1 | 1879 | { |
5443e65e | 1880 | // Parametrised "expected" error of the cluster reconstruction in Z |
1881 | ||
a9814c09 | 1882 | Double_t s = 9 * 9 /12.; |
5443e65e | 1883 | return s; |
1884 | } | |
1885 | ||
1886 | ||
1887 | //_____________________________________________________________________ | |
1888 | Double_t AliTRDtracker::GetX(Int_t sector, Int_t plane, Int_t local_tb) const | |
1889 | { | |
1890 | // | |
1891 | // Returns radial position which corresponds to time bin <local_tb> | |
1892 | // in tracking sector <sector> and plane <plane> | |
1893 | // | |
1894 | ||
1895 | Int_t index = fTrSec[sector]->CookTimeBinIndex(plane, local_tb); | |
1896 | Int_t pl = fTrSec[sector]->GetLayerNumber(index); | |
1897 | return fTrSec[sector]->GetLayer(pl)->GetX(); | |
1898 | ||
1899 | } | |
1900 | ||
1901 | ||
1902 | //_______________________________________________________ | |
1903 | AliTRDtracker::AliTRDpropagationLayer::AliTRDpropagationLayer(Double_t x, | |
a9814c09 | 1904 | Double_t dx, Double_t rho, Double_t rad_length, Int_t tb_index) |
5443e65e | 1905 | { |
0a29d0f1 | 1906 | // |
5443e65e | 1907 | // AliTRDpropagationLayer constructor |
0a29d0f1 | 1908 | // |
46d29e70 | 1909 | |
b3a5a838 | 1910 | fN = 0; fX = x; fdX = dx; fRho = rho; fX0 = rad_length; |
5443e65e | 1911 | fClusters = NULL; fIndex = NULL; fTimeBinIndex = tb_index; |
1912 | ||
46d29e70 | 1913 | |
5443e65e | 1914 | for(Int_t i=0; i < (Int_t) kZONES; i++) { |
1915 | fZc[i]=0; fZmax[i] = 0; | |
a819a5f7 | 1916 | } |
5443e65e | 1917 | |
1918 | fYmax = 0; | |
1919 | ||
1920 | if(fTimeBinIndex >= 0) { | |
ca6c93d6 | 1921 | fClusters = new AliTRDcluster*[kMAX_CLUSTER_PER_TIME_BIN]; |
1922 | fIndex = new UInt_t[kMAX_CLUSTER_PER_TIME_BIN]; | |
a819a5f7 | 1923 | } |
46d29e70 | 1924 | |
5443e65e | 1925 | fHole = kFALSE; |
1926 | fHoleZc = 0; | |
1927 | fHoleZmax = 0; | |
1928 | fHoleYc = 0; | |
1929 | fHoleYmax = 0; | |
1930 | fHoleRho = 0; | |
1931 | fHoleX0 = 0; | |
1932 | ||
1933 | } | |
1934 | ||
1935 | //_______________________________________________________ | |
1936 | void AliTRDtracker::AliTRDpropagationLayer::SetHole( | |
a9814c09 | 1937 | Double_t Zmax, Double_t Ymax, Double_t rho, |
1938 | Double_t rad_length, Double_t Yc, Double_t Zc) | |
5443e65e | 1939 | { |
1940 | // | |
1941 | // Sets hole in the layer | |
1942 | // | |
1943 | ||
1944 | fHole = kTRUE; | |
1945 | fHoleZc = Zc; | |
1946 | fHoleZmax = Zmax; | |
1947 | fHoleYc = Yc; | |
1948 | fHoleYmax = Ymax; | |
1949 | fHoleRho = rho; | |
b3a5a838 | 1950 | fHoleX0 = rad_length; |
5443e65e | 1951 | } |
1952 | ||
46d29e70 | 1953 | |
5443e65e | 1954 | //_______________________________________________________ |
1955 | AliTRDtracker::AliTRDtrackingSector::AliTRDtrackingSector(AliTRDgeometry* geo, Int_t gs, AliTRDparameter* par) | |
1956 | { | |
1957 | // | |
1958 | // AliTRDtrackingSector Constructor | |
1959 | // | |
1960 | ||
1961 | fGeom = geo; | |
1962 | fPar = par; | |
1963 | fGeomSector = gs; | |
1964 | fTzeroShift = 0.13; | |
1965 | fN = 0; | |
1966 | ||
1967 | for(UInt_t i=0; i < kMAX_TIME_BIN_INDEX; i++) fTimeBinIndex[i] = -1; | |
1968 | ||
1969 | ||
1970 | AliTRDpropagationLayer* ppl; | |
1971 | ||
b3a5a838 | 1972 | Double_t x, xin, xout, dx, rho, rad_length; |
5443e65e | 1973 | Int_t steps; |
46d29e70 | 1974 | |
5443e65e | 1975 | // set time bins in the gas of the TPC |
46d29e70 | 1976 | |
5443e65e | 1977 | xin = 246.055; xout = 254.055; steps = 20; dx = (xout-xin)/steps; |
b3a5a838 | 1978 | rho = 0.9e-3; rad_length = 28.94; |
5443e65e | 1979 | |
1980 | for(Int_t i=0; i<steps; i++) { | |
1981 | x = xin + i*dx + dx/2; | |
b3a5a838 | 1982 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,-1); |
5443e65e | 1983 | InsertLayer(ppl); |
46d29e70 | 1984 | } |
1985 | ||
5443e65e | 1986 | // set time bins in the outer field cage vessel |
46d29e70 | 1987 | |
b3a5a838 | 1988 | dx = 50e-4; xin = xout; xout = xin + dx; rho = 1.71; rad_length = 44.77; // Tedlar |
1989 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 1990 | InsertLayer(ppl); |
46d29e70 | 1991 | |
b3a5a838 | 1992 | dx = 0.02; xin = xout; xout = xin + dx; rho = 1.45; rad_length = 44.86; // prepreg |
1993 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 1994 | InsertLayer(ppl); |
1995 | ||
b3a5a838 | 1996 | dx = 2.; xin = xout; xout = xin + dx; rho = 1.45*0.02; rad_length = 41.28; // Nomex |
5443e65e | 1997 | steps = 5; dx = (xout - xin)/steps; |
1998 | for(Int_t i=0; i<steps; i++) { | |
1999 | x = xin + i*dx + dx/2; | |
b3a5a838 | 2000 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,-1); |
5443e65e | 2001 | InsertLayer(ppl); |
2002 | } | |
2003 | ||
b3a5a838 | 2004 | dx = 0.02; xin = xout; xout = xin + dx; rho = 1.45; rad_length = 44.86; // prepreg |
2005 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 2006 | InsertLayer(ppl); |
2007 | ||
b3a5a838 | 2008 | dx = 50e-4; xin = xout; xout = xin + dx; rho = 1.71; rad_length = 44.77; // Tedlar |
2009 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 2010 | InsertLayer(ppl); |
0a29d0f1 | 2011 | |
5443e65e | 2012 | |
2013 | // set time bins in CO2 | |
2014 | ||
2015 | xin = xout; xout = 275.0; | |
2016 | steps = 50; dx = (xout - xin)/steps; | |
b3a5a838 | 2017 | rho = 1.977e-3; rad_length = 36.2; |
5443e65e | 2018 | |
2019 | for(Int_t i=0; i<steps; i++) { | |
2020 | x = xin + i*dx + dx/2; | |
b3a5a838 | 2021 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,-1); |
5443e65e | 2022 | InsertLayer(ppl); |
2023 | } | |
2024 | ||
2025 | // set time bins in the outer containment vessel | |
2026 | ||
b3a5a838 | 2027 | dx = 50e-4; xin = xout; xout = xin + dx; rho = 2.7; rad_length = 24.01; // Al |
2028 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 2029 | InsertLayer(ppl); |
2030 | ||
b3a5a838 | 2031 | dx = 50e-4; xin = xout; xout = xin + dx; rho = 1.71; rad_length = 44.77; // Tedlar |
2032 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 2033 | InsertLayer(ppl); |
2034 | ||
b3a5a838 | 2035 | dx = 0.06; xin = xout; xout = xin + dx; rho = 1.45; rad_length = 44.86; // prepreg |
2036 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 2037 | InsertLayer(ppl); |
2038 | ||
b3a5a838 | 2039 | dx = 3.; xin = xout; xout = xin + dx; rho = 1.45*0.02; rad_length = 41.28; // Nomex |
5443e65e | 2040 | steps = 10; dx = (xout - xin)/steps; |
2041 | for(Int_t i=0; i<steps; i++) { | |
2042 | x = xin + i*dx + dx/2; | |
b3a5a838 | 2043 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,-1); |
5443e65e | 2044 | InsertLayer(ppl); |
2045 | } | |
2046 | ||
b3a5a838 | 2047 | dx = 0.06; xin = xout; xout = xin + dx; rho = 1.45; rad_length = 44.86; // prepreg |
2048 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 2049 | InsertLayer(ppl); |
2050 | ||
b3a5a838 | 2051 | dx = 50e-4; xin = xout; xout = xin + dx; rho = 1.71; rad_length = 44.77; // Tedlar |
2052 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 2053 | InsertLayer(ppl); |
2054 | ||
b3a5a838 | 2055 | dx = 50e-4; xin = xout; xout = xin + dx; rho = 2.7; rad_length = 24.01; // Al |
2056 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 2057 | InsertLayer(ppl); |
2058 | ||
2059 | Double_t xtrd = (Double_t) fGeom->Rmin(); | |
2060 | ||
2061 | // add layers between TPC and TRD (Air temporarily) | |
2062 | xin = xout; xout = xtrd; | |
2063 | steps = 50; dx = (xout - xin)/steps; | |
b3a5a838 | 2064 | rho = 1.2e-3; rad_length = 36.66; |
5443e65e | 2065 | |
2066 | for(Int_t i=0; i<steps; i++) { | |
2067 | x = xin + i*dx + dx/2; | |
b3a5a838 | 2068 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,-1); |
5443e65e | 2069 | InsertLayer(ppl); |
2070 | } | |
2071 | ||
2072 | ||
2073 | Double_t alpha=AliTRDgeometry::GetAlpha(); | |
2074 | ||
2075 | // add layers for each of the planes | |
2076 | ||
2077 | Double_t dxRo = (Double_t) fGeom->CroHght(); // Rohacell | |
2078 | Double_t dxSpace = (Double_t) fGeom->Cspace(); // Spacing between planes | |
2079 | Double_t dxAmp = (Double_t) fGeom->CamHght(); // Amplification region | |
2080 | Double_t dxDrift = (Double_t) fGeom->CdrHght(); // Drift region | |
2081 | Double_t dxRad = (Double_t) fGeom->CraHght(); // Radiator | |
2082 | Double_t dxTEC = dxRad + dxDrift + dxAmp + dxRo; | |
2083 | Double_t dxPlane = dxTEC + dxSpace; | |
2084 | ||
5443e65e | 2085 | Int_t tb, tb_index; |
2086 | const Int_t nChambers = AliTRDgeometry::Ncham(); | |
ca6c93d6 | 2087 | Double_t Ymax = 0, holeYmax = 0; |
2088 | Double_t * Zc = new Double_t[nChambers]; | |
2089 | Double_t * Zmax = new Double_t[nChambers]; | |
5443e65e | 2090 | Double_t holeZmax = 1000.; // the whole sector is missing |
2091 | ||
5443e65e | 2092 | for(Int_t plane = 0; plane < AliTRDgeometry::Nplan(); plane++) { |
2093 | ||
2094 | // Radiator | |
2095 | xin = xtrd + plane * dxPlane; xout = xin + dxRad; | |
b3a5a838 | 2096 | steps = 12; dx = (xout - xin)/steps; rho = 0.074; rad_length = 40.6; |
5443e65e | 2097 | for(Int_t i=0; i<steps; i++) { |
2098 | x = xin + i*dx + dx/2; | |
b3a5a838 | 2099 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,-1); |
5443e65e | 2100 | if((fGeom->GetPHOShole()) && (fGeomSector >= 2) && (fGeomSector <= 6)) { |
a9814c09 | 2101 | holeYmax = x*TMath::Tan(0.5*alpha); |
2102 | ppl->SetHole(holeYmax, holeZmax); | |
5443e65e | 2103 | } |
2104 | if((fGeom->GetRICHhole()) && (fGeomSector >= 12) && (fGeomSector <= 14)) { | |
a9814c09 | 2105 | holeYmax = x*TMath::Tan(0.5*alpha); |
2106 | ppl->SetHole(holeYmax, holeZmax); | |
5443e65e | 2107 | } |
2108 | InsertLayer(ppl); | |
2109 | } | |
2110 | ||
2111 | Ymax = fGeom->GetChamberWidth(plane)/2; | |
2112 | for(Int_t ch = 0; ch < nChambers; ch++) { | |
2113 | Zmax[ch] = fGeom->GetChamberLength(plane,ch)/2; | |
2114 | Float_t pad = fPar->GetRowPadSize(plane,ch,0); | |
2115 | Float_t row0 = fPar->GetRow0(plane,ch,0); | |
2116 | Int_t nPads = fPar->GetRowMax(plane,ch,0); | |
2117 | Zc[ch] = (pad * nPads)/2 + row0 - pad/2; | |
2118 | } | |
2119 | ||
2120 | dx = fPar->GetTimeBinSize(); | |
b3a5a838 | 2121 | rho = 0.00295 * 0.85; rad_length = 11.0; |
5443e65e | 2122 | |
2123 | Double_t x0 = (Double_t) fPar->GetTime0(plane); | |
2124 | Double_t xbottom = x0 - dxDrift; | |
2125 | Double_t xtop = x0 + dxAmp; | |
2126 | ||
2127 | // Amplification region | |
2128 | ||
2129 | steps = (Int_t) (dxAmp/dx); | |
2130 | ||
2131 | for(tb = 0; tb < steps; tb++) { | |
2132 | x = x0 + tb * dx + dx/2; | |
2133 | tb_index = CookTimeBinIndex(plane, -tb-1); | |
b3a5a838 | 2134 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,tb_index); |
5443e65e | 2135 | ppl->SetYmax(Ymax); |
2136 | for(Int_t ch = 0; ch < nChambers; ch++) { | |
a9814c09 | 2137 | ppl->SetZmax(ch, Zc[ch], Zmax[ch]); |
5443e65e | 2138 | } |
2139 | if((fGeom->GetPHOShole()) && (fGeomSector >= 2) && (fGeomSector <= 6)) { | |
a9814c09 | 2140 | holeYmax = x*TMath::Tan(0.5*alpha); |
2141 | ppl->SetHole(holeYmax, holeZmax); | |
5443e65e | 2142 | } |
2143 | if((fGeom->GetRICHhole()) && (fGeomSector >= 12) && (fGeomSector <= 14)) { | |
a9814c09 | 2144 | holeYmax = x*TMath::Tan(0.5*alpha); |
2145 | ppl->SetHole(holeYmax, holeZmax); | |
5443e65e | 2146 | } |
2147 | InsertLayer(ppl); | |
2148 | } | |
2149 | tb_index = CookTimeBinIndex(plane, -steps); | |
2150 | x = (x + dx/2 + xtop)/2; | |
2151 | dx = 2*(xtop-x); | |
b3a5a838 | 2152 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,tb_index); |
5443e65e | 2153 | ppl->SetYmax(Ymax); |
2154 | for(Int_t ch = 0; ch < nChambers; ch++) { | |
2155 | ppl->SetZmax(ch, Zc[ch], Zmax[ch]); | |
2156 | } | |
2157 | if((fGeom->GetPHOShole()) && (fGeomSector >= 2) && (fGeomSector <= 6)) { | |
2158 | holeYmax = x*TMath::Tan(0.5*alpha); | |
2159 | ppl->SetHole(holeYmax, holeZmax); | |
2160 | } | |
2161 | if((fGeom->GetRICHhole()) && (fGeomSector >= 12) && (fGeomSector <= 14)) { | |
2162 | holeYmax = x*TMath::Tan(0.5*alpha); | |
2163 | ppl->SetHole(holeYmax, holeZmax); | |
2164 | } | |
2165 | InsertLayer(ppl); | |
2166 | ||
2167 | // Drift region | |
2168 | dx = fPar->GetTimeBinSize(); | |
2169 | steps = (Int_t) (dxDrift/dx); | |
2170 | ||
2171 | for(tb = 0; tb < steps; tb++) { | |
2172 | x = x0 - tb * dx - dx/2; | |
2173 | tb_index = CookTimeBinIndex(plane, tb); | |
2174 | ||
b3a5a838 | 2175 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,tb_index); |
5443e65e | 2176 | ppl->SetYmax(Ymax); |
2177 | for(Int_t ch = 0; ch < nChambers; ch++) { | |
a9814c09 | 2178 | ppl->SetZmax(ch, Zc[ch], Zmax[ch]); |
5443e65e | 2179 | } |
2180 | if((fGeom->GetPHOShole()) && (fGeomSector >= 2) && (fGeomSector <= 6)) { | |
a9814c09 | 2181 | holeYmax = x*TMath::Tan(0.5*alpha); |
2182 | ppl->SetHole(holeYmax, holeZmax); | |
5443e65e | 2183 | } |
2184 | if((fGeom->GetRICHhole()) && (fGeomSector >= 12) && (fGeomSector <= 14)) { | |
a9814c09 | 2185 | holeYmax = x*TMath::Tan(0.5*alpha); |
2186 | ppl->SetHole(holeYmax, holeZmax); | |
5443e65e | 2187 | } |
2188 | InsertLayer(ppl); | |
2189 | } | |
2190 | tb_index = CookTimeBinIndex(plane, steps); | |
2191 | x = (x - dx/2 + xbottom)/2; | |
2192 | dx = 2*(x-xbottom); | |
b3a5a838 | 2193 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,tb_index); |
5443e65e | 2194 | ppl->SetYmax(Ymax); |
2195 | for(Int_t ch = 0; ch < nChambers; ch++) { | |
2196 | ppl->SetZmax(ch, Zc[ch], Zmax[ch]); | |
2197 | } | |
2198 | if((fGeom->GetPHOShole()) && (fGeomSector >= 2) && (fGeomSector <= 6)) { | |
2199 | holeYmax = x*TMath::Tan(0.5*alpha); | |
2200 | ppl->SetHole(holeYmax, holeZmax); | |
2201 | } | |
2202 | if((fGeom->GetRICHhole()) && (fGeomSector >= 12) && (fGeomSector <= 14)) { | |
2203 | holeYmax = x*TMath::Tan(0.5*alpha); | |
2204 | ppl->SetHole(holeYmax, holeZmax); | |
2205 | } | |
2206 | InsertLayer(ppl); | |
2207 | ||
2208 | // Pad Plane | |
b3a5a838 | 2209 | xin = xtop; dx = 0.025; xout = xin + dx; rho = 1.7; rad_length = 33.0; |
2210 | ppl = new AliTRDpropagationLayer(xin+dx/2,dx,rho,rad_length,-1); | |
5443e65e | 2211 | if((fGeom->GetPHOShole()) && (fGeomSector >= 2) && (fGeomSector <= 6)) { |
2212 | holeYmax = (xin+dx/2)*TMath::Tan(0.5*alpha); | |
2213 | ppl->SetHole(holeYmax, holeZmax); | |
2214 | } | |
2215 | if((fGeom->GetRICHhole()) && (fGeomSector >= 12) && (fGeomSector <= 14)) { | |
2216 | holeYmax = (xin+dx/2)*TMath::Tan(0.5*alpha); | |
2217 | ppl->SetHole(holeYmax, holeZmax); | |
2218 | } | |
2219 | InsertLayer(ppl); | |
2220 | ||
2221 | // Rohacell | |
2222 | xin = xout; xout = xtrd + (plane + 1) * dxPlane - dxSpace; | |
b3a5a838 | 2223 | steps = 5; dx = (xout - xin)/steps; rho = 0.074; rad_length = 40.6; |
5443e65e | 2224 | for(Int_t i=0; i<steps; i++) { |
2225 | x = xin + i*dx + dx/2; | |
b3a5a838 | 2226 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,-1); |
5443e65e | 2227 | if((fGeom->GetPHOShole()) && (fGeomSector >= 2) && (fGeomSector <= 6)) { |
a9814c09 | 2228 | holeYmax = x*TMath::Tan(0.5*alpha); |
2229 | ppl->SetHole(holeYmax, holeZmax); | |
5443e65e | 2230 | } |
2231 | if((fGeom->GetRICHhole()) && (fGeomSector >= 12) && (fGeomSector <= 14)) { | |
a9814c09 | 2232 | holeYmax = x*TMath::Tan(0.5*alpha); |
2233 | ppl->SetHole(holeYmax, holeZmax); | |
5443e65e | 2234 | } |
2235 | InsertLayer(ppl); | |
2236 | } | |
2237 | ||
2238 | // Space between the chambers, air | |
2239 | xin = xout; xout = xtrd + (plane + 1) * dxPlane; | |
b3a5a838 | 2240 | steps = 5; dx = (xout - xin)/steps; rho = 1.29e-3; rad_length = 36.66; |
5443e65e | 2241 | for(Int_t i=0; i<steps; i++) { |
2242 | x = xin + i*dx + dx/2; | |
b3a5a838 | 2243 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,-1); |
5443e65e | 2244 | if((fGeom->GetPHOShole()) && (fGeomSector >= 2) && (fGeomSector <= 6)) { |
a9814c09 | 2245 | holeYmax = x*TMath::Tan(0.5*alpha); |
2246 | ppl->SetHole(holeYmax, holeZmax); | |
5443e65e | 2247 | } |
2248 | if((fGeom->GetRICHhole()) && (fGeomSector >= 12) && (fGeomSector <= 14)) { | |
a9814c09 | 2249 | holeYmax = x*TMath::Tan(0.5*alpha); |
2250 | ppl->SetHole(holeYmax, holeZmax); | |
5443e65e | 2251 | } |
2252 | InsertLayer(ppl); | |
2253 | } | |
2254 | } | |
2255 | ||
2256 | // Space between the TRD and RICH | |
2257 | Double_t xRICH = 500.; | |
2258 | xin = xout; xout = xRICH; | |
b3a5a838 | 2259 | steps = 200; dx = (xout - xin)/steps; rho = 1.29e-3; rad_length = 36.66; |
5443e65e | 2260 | for(Int_t i=0; i<steps; i++) { |
2261 | x = xin + i*dx + dx/2; | |
b3a5a838 | 2262 | ppl = new AliTRDpropagationLayer(x,dx,rho,rad_length,-1); |
5443e65e | 2263 | InsertLayer(ppl); |
2264 | } | |
2265 | ||
2266 | MapTimeBinLayers(); | |
ca6c93d6 | 2267 | delete [] Zc; |
2268 | delete [] Zmax; | |
5443e65e | 2269 | |
2270 | } | |
2271 | ||
2272 | //______________________________________________________ | |
2273 | ||
2274 | Int_t AliTRDtracker::AliTRDtrackingSector::CookTimeBinIndex(Int_t plane, Int_t local_tb) const | |
2275 | { | |
2276 | // | |
2277 | // depending on the digitization parameters calculates "global" | |
2278 | // time bin index for timebin <local_tb> in plane <plane> | |
2279 | // | |
2280 | ||
2281 | Double_t dxAmp = (Double_t) fGeom->CamHght(); // Amplification region | |
2282 | Double_t dxDrift = (Double_t) fGeom->CdrHght(); // Drift region | |
2283 | Double_t dx = (Double_t) fPar->GetTimeBinSize(); | |
2284 | ||
2285 | Int_t tbAmp = fPar->GetTimeBefore(); | |
2286 | Int_t maxAmp = (Int_t) ((dxAmp+0.000001)/dx); | |
b3a5a838 | 2287 | if(kTRUE) maxAmp = 0; // intentional until we change parameter class |
5443e65e | 2288 | Int_t tbDrift = fPar->GetTimeMax(); |
2289 | Int_t maxDrift = (Int_t) ((dxDrift+0.000001)/dx); | |
2290 | ||
2291 | Int_t tb_per_plane = TMath::Min(tbAmp,maxAmp) + TMath::Min(tbDrift,maxDrift); | |
2292 | ||
b3a5a838 | 2293 | Int_t gtb = (plane+1) * tb_per_plane - local_tb - 1 - TMath::Min(tbAmp,maxAmp); |
5443e65e | 2294 | |
2295 | if((local_tb < 0) && | |
2296 | (TMath::Abs(local_tb) > TMath::Min(tbAmp,maxAmp))) return -1; | |
2297 | if(local_tb >= TMath::Min(tbDrift,maxDrift)) return -1; | |
2298 | ||
2299 | return gtb; | |
2300 | ||
2301 | ||
2302 | } | |
2303 | ||
2304 | //______________________________________________________ | |
2305 | ||
2306 | void AliTRDtracker::AliTRDtrackingSector::MapTimeBinLayers() | |
2307 | { | |
2308 | // | |
2309 | // For all sensitive time bins sets corresponding layer index | |
2310 | // in the array fTimeBins | |
2311 | // | |
2312 | ||
2313 | Int_t index; | |
2314 | ||
2315 | for(Int_t i = 0; i < fN; i++) { | |
2316 | index = fLayers[i]->GetTimeBinIndex(); | |
2317 | ||
2318 | // printf("gtb %d -> pl %d -> x %f \n", index, i, fLayers[i]->GetX()); | |
2319 | ||
2320 | if(index < 0) continue; | |
2321 | if(index >= (Int_t) kMAX_TIME_BIN_INDEX) { | |
2322 | printf("*** AliTRDtracker::MapTimeBinLayers: \n"); | |
2323 | printf(" index %d exceeds allowed maximum of %d!\n", | |
a9814c09 | 2324 | index, kMAX_TIME_BIN_INDEX-1); |
5443e65e | 2325 | continue; |
2326 | } | |
2327 | fTimeBinIndex[index] = i; | |
2328 | } | |
2329 | ||
2330 | Double_t x1, dx1, x2, dx2, gap; | |
2331 | ||
2332 | for(Int_t i = 0; i < fN-1; i++) { | |
2333 | x1 = fLayers[i]->GetX(); | |
2334 | dx1 = fLayers[i]->GetdX(); | |
2335 | x2 = fLayers[i+1]->GetX(); | |
2336 | dx2 = fLayers[i+1]->GetdX(); | |
2337 | gap = (x2 - dx2/2) - (x1 + dx1/2); | |
2338 | if(gap < -0.01) { | |
2339 | printf("*** warning: layers %d and %d are overlayed:\n",i,i+1); | |
2340 | printf(" %f + %f + %f > %f\n", x1, dx1/2, dx2/2, x2); | |
2341 | } | |
2342 | if(gap > 0.01) { | |
2343 | printf("*** warning: layers %d and %d have a large gap:\n",i,i+1); | |
2344 | printf(" (%f - %f) - (%f + %f) = %f\n", | |
a9814c09 | 2345 | x2, dx2/2, x1, dx1, gap); |
5443e65e | 2346 | } |
2347 | } | |
2348 | } | |
2349 | ||
2350 | ||
2351 | //______________________________________________________ | |
2352 | ||
2353 | ||
2354 | Int_t AliTRDtracker::AliTRDtrackingSector::GetLayerNumber(Double_t x) const | |
2355 | { | |
2356 | // | |
2357 | // Returns the number of time bin which in radial position is closest to <x> | |
2358 | // | |
2359 | ||
2360 | if(x >= fLayers[fN-1]->GetX()) return fN-1; | |
2361 | if(x <= fLayers[0]->GetX()) return 0; | |
2362 | ||
2363 | Int_t b=0, e=fN-1, m=(b+e)/2; | |
2364 | for (; b<e; m=(b+e)/2) { | |
2365 | if (x > fLayers[m]->GetX()) b=m+1; | |
2366 | else e=m; | |
2367 | } | |
2368 | if(TMath::Abs(x - fLayers[m]->GetX()) > | |
2369 | TMath::Abs(x - fLayers[m+1]->GetX())) return m+1; | |
2370 | else return m; | |
2371 | ||
2372 | } | |
2373 | ||
2374 | //______________________________________________________ | |
2375 | ||
2376 | Int_t AliTRDtracker::AliTRDtrackingSector::GetInnerTimeBin() const | |
2377 | { | |
2378 | // | |
2379 | // Returns number of the innermost SENSITIVE propagation layer | |
2380 | // | |
2381 | ||
2382 | return GetLayerNumber(0); | |
2383 | } | |
2384 | ||
2385 | //______________________________________________________ | |
2386 | ||
2387 | Int_t AliTRDtracker::AliTRDtrackingSector::GetOuterTimeBin() const | |
2388 | { | |
2389 | // | |
2390 | // Returns number of the outermost SENSITIVE time bin | |
2391 | // | |
2392 | ||
2393 | return GetLayerNumber(GetNumberOfTimeBins() - 1); | |
46d29e70 | 2394 | } |
2395 | ||
5443e65e | 2396 | //______________________________________________________ |
2397 | ||
2398 | Int_t AliTRDtracker::AliTRDtrackingSector::GetNumberOfTimeBins() const | |
2399 | { | |
2400 | // | |
2401 | // Returns number of SENSITIVE time bins | |
2402 | // | |
2403 | ||
2404 | Int_t tb, layer; | |
2405 | for(tb = kMAX_TIME_BIN_INDEX-1; tb >=0; tb--) { | |
2406 | layer = GetLayerNumber(tb); | |
2407 | if(layer>=0) break; | |
2408 | } | |
2409 | return tb+1; | |
2410 | } | |
2411 | ||
2412 | //______________________________________________________ | |
2413 | ||
2414 | void AliTRDtracker::AliTRDtrackingSector::InsertLayer(AliTRDpropagationLayer* pl) | |
2415 | { | |
2416 | // | |
2417 | // Insert layer <pl> in fLayers array. | |
2418 | // Layers are sorted according to X coordinate. | |
2419 | ||
2420 | if ( fN == ((Int_t) kMAX_LAYERS_PER_SECTOR)) { | |
2421 | printf("AliTRDtrackingSector::InsertLayer(): Too many layers !\n"); | |
2422 | return; | |
2423 | } | |
2424 | if (fN==0) {fLayers[fN++] = pl; return;} | |
2425 | Int_t i=Find(pl->GetX()); | |
2426 | ||
2427 | memmove(fLayers+i+1 ,fLayers+i,(fN-i)*sizeof(AliTRDpropagationLayer*)); | |
2428 | fLayers[i]=pl; fN++; | |
2429 | ||
2430 | } | |
2431 | ||
2432 | //______________________________________________________ | |
2433 | ||
2434 | Int_t AliTRDtracker::AliTRDtrackingSector::Find(Double_t x) const | |
2435 | { | |
2436 | // | |
2437 | // Returns index of the propagation layer nearest to X | |
2438 | // | |
2439 | ||
2440 | if (x <= fLayers[0]->GetX()) return 0; | |
2441 | if (x > fLayers[fN-1]->GetX()) return fN; | |
2442 | Int_t b=0, e=fN-1, m=(b+e)/2; | |
2443 | for (; b<e; m=(b+e)/2) { | |
2444 | if (x > fLayers[m]->GetX()) b=m+1; | |
2445 | else e=m; | |
2446 | } | |
2447 | return m; | |
2448 | } | |
2449 | ||
2450 | //______________________________________________________ | |
2451 | ||
2452 | void AliTRDtracker::AliTRDpropagationLayer::GetPropagationParameters( | |
b3a5a838 | 2453 | Double_t y, Double_t z, Double_t &dx, Double_t &rho, Double_t &rad_length, |
a9814c09 | 2454 | Bool_t &lookForCluster) const |
5443e65e | 2455 | { |
2456 | // | |
b3a5a838 | 2457 | // Returns radial step <dx>, density <rho>, rad. length <rad_length>, |
5443e65e | 2458 | // and sensitivity <lookForCluster> in point <y,z> |
2459 | // | |
2460 | ||
2461 | dx = fdX; | |
2462 | rho = fRho; | |
b3a5a838 | 2463 | rad_length = fX0; |
5443e65e | 2464 | lookForCluster = kFALSE; |
2465 | ||
2466 | // check dead regions | |
2467 | if(fTimeBinIndex >= 0) { | |
2468 | for(Int_t ch = 0; ch < (Int_t) kZONES; ch++) { | |
2469 | if(TMath::Abs(z - fZc[ch]) < fZmax[ch]) | |
a9814c09 | 2470 | lookForCluster = kTRUE; |
2471 | // else { rho = 1.7; rad_length = 33.0; } // G10 | |
5443e65e | 2472 | } |
2473 | if(TMath::Abs(y) > fYmax) lookForCluster = kFALSE; | |
2474 | if(!lookForCluster) { | |
b3a5a838 | 2475 | // rho = 1.7; rad_length = 33.0; // G10 |
5443e65e | 2476 | } |
2477 | } | |
2478 | ||
2479 | // check hole | |
2480 | if(fHole && (TMath::Abs(y - fHoleYc) < fHoleYmax) && | |
2481 | (TMath::Abs(z - fHoleZc) < fHoleZmax)) { | |
2482 | lookForCluster = kFALSE; | |
2483 | rho = fHoleRho; | |
b3a5a838 | 2484 | rad_length = fHoleX0; |
5443e65e | 2485 | } |
2486 | ||
2487 | return; | |
2488 | } | |
2489 | ||
2490 | //______________________________________________________ | |
2491 | ||
2492 | void AliTRDtracker::AliTRDpropagationLayer::InsertCluster(AliTRDcluster* c, | |
a9814c09 | 2493 | UInt_t index) { |
5443e65e | 2494 | |
2495 | // Insert cluster in cluster array. | |
2496 | // Clusters are sorted according to Y coordinate. | |
2497 | ||
2498 | if(fTimeBinIndex < 0) { | |
2499 | printf("*** attempt to insert cluster into non-sensitive time bin!\n"); | |
2500 | return; | |
2501 | } | |
2502 | ||
2503 | if (fN== (Int_t) kMAX_CLUSTER_PER_TIME_BIN) { | |
2504 | printf("AliTRDpropagationLayer::InsertCluster(): Too many clusters !\n"); | |
2505 | return; | |
2506 | } | |
2507 | if (fN==0) {fIndex[0]=index; fClusters[fN++]=c; return;} | |
2508 | Int_t i=Find(c->GetY()); | |
2509 | memmove(fClusters+i+1 ,fClusters+i,(fN-i)*sizeof(AliTRDcluster*)); | |
2510 | memmove(fIndex +i+1 ,fIndex +i,(fN-i)*sizeof(UInt_t)); | |
2511 | fIndex[i]=index; fClusters[i]=c; fN++; | |
2512 | } | |
2513 | ||
2514 | //______________________________________________________ | |
2515 | ||
2516 | Int_t AliTRDtracker::AliTRDpropagationLayer::Find(Double_t y) const { | |
2517 | ||
2518 | // Returns index of the cluster nearest in Y | |
2519 | ||
2520 | if (y <= fClusters[0]->GetY()) return 0; | |
2521 | if (y > fClusters[fN-1]->GetY()) return fN; | |
2522 | Int_t b=0, e=fN-1, m=(b+e)/2; | |
2523 | for (; b<e; m=(b+e)/2) { | |
2524 | if (y > fClusters[m]->GetY()) b=m+1; | |
2525 | else e=m; | |
2526 | } | |
2527 | return m; | |
2528 | } | |
2529 | ||
fd621f36 | 2530 | //--------------------------------------------------------- |
5443e65e | 2531 | |
fd621f36 | 2532 | Double_t AliTRDtracker::GetTiltFactor(const AliTRDcluster* c) { |
2533 | // | |
2534 | // Returns correction factor for tilted pads geometry | |
2535 | // | |
5443e65e | 2536 | |
fd621f36 | 2537 | Double_t h01 = sin(TMath::Pi() / 180.0 * fPar->GetTiltingAngle()); |
2538 | Int_t det = c->GetDetector(); | |
2539 | Int_t plane = fGeom->GetPlane(det); | |
b3a5a838 | 2540 | |
2541 | if((plane == 1) || (plane == 3) || (plane == 5)) h01=-h01; | |
b8dc2353 | 2542 | |
2543 | if(fNoTilt) h01 = 0; | |
fd621f36 | 2544 | |
2545 | return h01; | |
2546 | } | |
5443e65e | 2547 | |
2548 | ||
2549 | ||
2550 | ||
2551 | ||
2552 |