]>
Commit | Line | Data |
---|---|---|
b1886074 | 1 | //$Id$ |
2 | ||
3 | // Author: Anders Vestbo <mailto:vestbo@fi.uib.no> | |
4 | //*-- Copyright © ASV | |
d712b5b8 | 5 | |
6 | #include <math.h> | |
7 | ||
8 | #include "AliL3Logging.h" | |
9 | #include "AliL3Transform.h" | |
10 | #include "AliL3Defs.h" | |
11 | #include "AliL3HoughTrack.h" | |
12 | ||
b1886074 | 13 | //_____________________________________________________________ |
14 | // AliL3HoughTrack | |
15 | // | |
16 | // Track class for Hough tracklets | |
17 | ||
d712b5b8 | 18 | ClassImp(AliL3HoughTrack) |
19 | ||
20 | ||
21 | AliL3HoughTrack::AliL3HoughTrack() | |
22 | { | |
23 | //Constructor | |
24 | ||
25 | fWeight = 0; | |
26 | fMinDist=0; | |
27 | fTransform = new AliL3Transform(); | |
28 | fDLine = 0; | |
29 | fPsiLine = 0; | |
30 | fIsHelix = true; | |
31 | fEtaIndex = -1; | |
32 | } | |
33 | ||
34 | ||
35 | AliL3HoughTrack::~AliL3HoughTrack() | |
36 | { | |
37 | //Destructor | |
38 | if(fTransform) | |
39 | delete fTransform; | |
40 | } | |
41 | ||
42 | void AliL3HoughTrack::Set(AliL3Track *track) | |
43 | { | |
44 | ||
45 | AliL3HoughTrack *tpt = (AliL3HoughTrack*)track; | |
46 | SetTrackParameters(tpt->GetKappa(),tpt->GetPhi0(),tpt->GetWeight()); | |
47 | SetEtaIndex(tpt->GetEtaIndex()); | |
48 | SetEta(tpt->GetEta()); | |
ae4d5136 | 49 | SetTgl(tpt->GetTgl()); |
d712b5b8 | 50 | SetPsi(tpt->GetPsi()); |
51 | SetCenterX(tpt->GetCenterX()); | |
52 | SetCenterY(tpt->GetCenterY()); | |
53 | SetFirstPoint(tpt->GetFirstPointX(),tpt->GetFirstPointY(),tpt->GetFirstPointZ()); | |
54 | SetLastPoint(tpt->GetLastPointX(),tpt->GetLastPointY(),tpt->GetLastPointZ()); | |
55 | SetCharge(tpt->GetCharge()); | |
56 | SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow()); | |
b1886074 | 57 | SetSlice(tpt->GetSlice()); |
242c143e | 58 | SetNHits(1); |
d712b5b8 | 59 | return; |
60 | ||
61 | fWeight = tpt->GetWeight(); | |
62 | fDLine = tpt->GetDLine(); | |
63 | fPsiLine = tpt->GetPsiLine(); | |
64 | SetNHits(tpt->GetWeight()); | |
65 | SetRowRange(tpt->GetFirstRow(),tpt->GetLastRow()); | |
66 | fIsHelix = false; | |
67 | ||
68 | ||
69 | } | |
70 | ||
71 | Int_t AliL3HoughTrack::Compare(const AliL3Track *tpt) const | |
72 | { | |
73 | AliL3HoughTrack *track = (AliL3HoughTrack*)tpt; | |
74 | if(track->GetWeight() < GetWeight()) return 1; | |
75 | if(track->GetWeight() > GetWeight()) return -1; | |
76 | return 0; | |
77 | } | |
78 | ||
79 | void AliL3HoughTrack::SetEta(Double_t f) | |
80 | { | |
b1886074 | 81 | //Set eta, and calculate fTanl, which is the tan of dipangle |
82 | ||
d712b5b8 | 83 | fEta = f; |
84 | Double_t theta = 2*atan(exp(-1.*fEta)); | |
8b007a2d | 85 | Double_t dipangle = Pi/2 - theta; |
86 | Double_t tgl = tan(dipangle); | |
d712b5b8 | 87 | SetTgl(tgl); |
88 | } | |
89 | ||
90 | void AliL3HoughTrack::CalculateHelix() | |
91 | { | |
92 | return; | |
93 | } | |
94 | ||
95 | void AliL3HoughTrack::UpdateToFirstRow() | |
96 | { | |
97 | //Update the track parameters to the point where track cross | |
98 | //its first padrow. | |
99 | ||
100 | //Get the crossing point with the first padrow: | |
101 | Float_t xyz[3]; | |
102 | if(!GetCrossingPoint(GetFirstRow(),xyz)) | |
103 | LOG(AliL3Log::kWarning,"AliL3HoughTrack::UpdateToFirstRow()","Track parameters") | |
104 | <<AliL3Log::kDec<<"Track does not cross padrow "<<GetFirstRow()<<" centerx " | |
105 | <<GetCenterX()<<" centery "<<GetCenterY()<<" Radius "<<GetRadius()<<" tgl "<<GetTgl()<<ENDLOG; | |
106 | ||
8b007a2d | 107 | //printf("Track with eta %f tgl %f crosses at x %f y %f z %f on padrow %d\n",GetEta(),GetTgl(),xyz[0],xyz[1],xyz[2],GetFirstRow()); |
108 | //printf("Before: first %f %f %f tgl %f center %f %f charge %d\n",GetFirstPointX(),GetFirstPointY(),GetFirstPointZ(),GetTgl(),GetCenterX(),GetCenterY(),GetCharge()); | |
109 | ||
d712b5b8 | 110 | Double_t radius = sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]); |
111 | ||
112 | //Get the track parameters | |
b1886074 | 113 | |
114 | /* | |
115 | Double_t x0 = GetR0() * cos(GetPhi0()) ; | |
116 | Double_t y0 = GetR0() * sin(GetPhi0()) ; | |
117 | */ | |
118 | Double_t rc = GetRadius();//fabs(GetPt()) / ( BFACT * BField ) ; | |
d712b5b8 | 119 | Double_t tPhi0 = GetPsi() + GetCharge() * 0.5 * pi / fabs(GetCharge()) ; |
b1886074 | 120 | Double_t xc = GetCenterX();//x0 - rc * cos(tPhi0) ; |
121 | Double_t yc = GetCenterY();//y0 - rc * sin(tPhi0) ; | |
d712b5b8 | 122 | |
123 | //Check helix and cylinder intersect | |
124 | Double_t fac1 = xc*xc + yc*yc ; | |
125 | Double_t sfac = sqrt( fac1 ) ; | |
b1886074 | 126 | |
d712b5b8 | 127 | if ( fabs(sfac-rc) > radius || fabs(sfac+rc) < radius ) { |
128 | LOG(AliL3Log::kError,"AliL3HoughTrack::UpdateToFirstRow","Tracks")<<AliL3Log::kDec<< | |
129 | "Track does not intersect"<<ENDLOG; | |
130 | return; | |
131 | } | |
132 | ||
133 | //Find intersection | |
134 | Double_t fac2 = ( radius*radius + fac1 - rc*rc) / (2.00 * radius * sfac ) ; | |
135 | Double_t phi = atan2(yc,xc) + GetCharge()*acos(fac2) ; | |
136 | Double_t td = atan2(radius*sin(phi) - yc,radius*cos(phi) - xc) ; | |
137 | ||
138 | //Intersection in z | |
139 | if ( td < 0 ) td = td + 2. * pi ; | |
140 | Double_t deltat = fmod((-GetCharge()*td + GetCharge()*tPhi0),2*pi) ; | |
141 | if ( deltat < 0. ) deltat += 2. * pi ; | |
142 | if ( deltat > 2.*pi ) deltat -= 2. * pi ; | |
143 | Double_t z = GetZ0() + rc * GetTgl() * deltat ; | |
8b007a2d | 144 | |
d712b5b8 | 145 | |
146 | Double_t xExtra = radius * cos(phi) ; | |
147 | Double_t yExtra = radius * sin(phi) ; | |
148 | ||
149 | Double_t tPhi = atan2(yExtra-yc,xExtra-xc); | |
150 | ||
151 | //if ( tPhi < 0 ) tPhi += 2. * M_PI ; | |
152 | ||
153 | Double_t tPsi = tPhi - GetCharge() * 0.5 * pi / fabs(GetCharge()) ; | |
154 | if ( tPsi > 2. * pi ) tPsi -= 2. * pi ; | |
155 | if ( tPsi < 0. ) tPsi += 2. * pi ; | |
156 | ||
157 | //And finally, update the track parameters | |
d712b5b8 | 158 | SetR0(radius); |
159 | SetPhi0(phi); | |
160 | SetZ0(z); | |
161 | SetPsi(tPsi); | |
162 | SetFirstPoint(xyz[0],xyz[1],z); | |
8b007a2d | 163 | //printf("After: first %f %f %f tgl %f center %f %f charge %d\n",GetFirstPointX(),GetFirstPointY(),GetFirstPointZ(),GetTgl(),GetCenterX(),GetCenterY(),GetCharge()); |
164 | ||
165 | //printf("First point set %f %f %f\n",xyz[0],xyz[1],z); | |
d712b5b8 | 166 | |
167 | //Also, set the coordinates of the point where track crosses last padrow: | |
168 | GetCrossingPoint(GetLastRow(),xyz); | |
169 | SetLastPoint(xyz[0],xyz[1],xyz[2]); | |
8b007a2d | 170 | //printf("last point %f %f %f\n",xyz[0],xyz[1],xyz[2]); |
d712b5b8 | 171 | } |
172 | ||
173 | void AliL3HoughTrack::SetTrackParameters(Double_t kappa,Double_t phi,Int_t weight) | |
174 | { | |
175 | ||
176 | fWeight = weight; | |
177 | fMinDist = 100000; | |
178 | SetKappa(kappa); | |
179 | SetPhi0(phi); | |
b1886074 | 180 | Double_t pt = fabs(BFACT*BField/kappa); |
d712b5b8 | 181 | SetPt(pt); |
182 | Double_t radius = 1/fabs(kappa); | |
183 | SetRadius(radius); | |
184 | SetFirstPoint(0,0,0); | |
8b007a2d | 185 | SetPsi(phi); //Psi = Phi when first point is vertex |
186 | SetR0(0); | |
d712b5b8 | 187 | Double_t charge = -1.*kappa; |
8b007a2d | 188 | SetCharge((Int_t)copysign(1.,charge)); |
189 | ||
190 | Double_t trackPhi0 = GetPsi() + charge*0.5*Pi/fabs(charge); | |
191 | Double_t xc = GetFirstPointX() - GetRadius() * cos(trackPhi0) ; | |
192 | Double_t yc = GetFirstPointY() - GetRadius() * sin(trackPhi0) ; | |
d712b5b8 | 193 | SetCenterX(xc); |
194 | SetCenterY(yc); | |
242c143e | 195 | SetNHits(1); //just for the trackarray IO |
d712b5b8 | 196 | fIsHelix = true; |
8b007a2d | 197 | |
242c143e | 198 | |
d712b5b8 | 199 | } |
200 | ||
201 | void AliL3HoughTrack::SetLineParameters(Double_t psi,Double_t D,Int_t weight,Int_t *rowrange,Int_t ref_row) | |
202 | { | |
203 | //Initialize a track piece, not yet a track | |
204 | //Used in case of straight line transformation | |
205 | ||
206 | //Transform line parameters to coordinate system of slice: | |
207 | ||
208 | D = D + fTransform->Row2X(ref_row)*cos(psi); | |
209 | ||
210 | fDLine = D; | |
211 | fPsiLine = psi; | |
212 | fWeight = weight; | |
213 | SetNHits(weight); | |
214 | SetRowRange(rowrange[0],rowrange[1]); | |
215 | fIsHelix = false; | |
216 | ||
217 | } | |
218 | ||
219 | void AliL3HoughTrack::SetBestMCid(Int_t mcid,Double_t min_dist) | |
220 | { | |
221 | ||
222 | if(min_dist < fMinDist) | |
223 | { | |
224 | fMinDist = min_dist; | |
225 | SetMCid(mcid); | |
226 | } | |
227 | ||
228 | } | |
229 | ||
230 | void AliL3HoughTrack::GetLineCrossingPoint(Int_t padrow,Double_t *xy) | |
231 | { | |
232 | ||
233 | ||
234 | if(fIsHelix) | |
235 | { | |
236 | printf("AliL3HoughTrack::GetLineCrossingPoint : Track is not a line\n"); | |
237 | return; | |
238 | } | |
239 | ||
240 | Double_t xhit = fTransform->Row2X(padrow); | |
241 | Double_t a = -1/tan(fPsiLine); | |
242 | Double_t b = fDLine/sin(fPsiLine); | |
243 | ||
244 | Double_t yhit = a*xhit + b; | |
245 | xy[0] = xhit; | |
246 | xy[1] = yhit; | |
247 | ||
248 | } | |
249 | ||
250 | /* | |
251 | Double_t AliL3HoughTrack::GetCrossingAngle(Int_t padrow) | |
252 | { | |
253 | //Calculate the crossing angle between track and given padrow. | |
254 | ||
255 | if(!fIsHelix) | |
256 | { | |
257 | printf("AliL3HoughTrack::GetCrossingAngle : Track is not a helix\n"); | |
258 | return 0; | |
259 | } | |
260 | ||
261 | if(!IsLocal()) | |
262 | { | |
263 | printf("Track is not given in local coordinates\n"); | |
264 | return 0; | |
265 | } | |
266 | ||
267 | Float_t xyz[3]; | |
268 | if(!GetCrossingPoint(padrow,xyz)) | |
269 | printf("AliL3HoughTrack::GetCrossingPoint : Track does not cross line!!\n"); | |
270 | ||
271 | ||
272 | //Convert center of curvature to local coordinates: | |
273 | //Float_t xyz_coc[3] = {GetCenterX(),GetCenterY(),0}; | |
274 | //fTransform->Global2Local(xyz_coc,slice); | |
275 | ||
276 | //Take the dot product of the tangent vector of the track, and | |
277 | //vector perpendicular to the padrow. | |
278 | ||
279 | Double_t tangent[2]; | |
280 | //tangent[1] = (xyz[0] - xyz_coc[0])/GetRadius(); | |
281 | //tangent[0] = -1.*(xyz[1] - xyz_coc[1])/GetRadius(); | |
282 | tangent[1] = (xyz[0] - GetCenterX())/GetRadius(); | |
283 | tangent[0] = -1.*(xyz[1] - GetCenterY())/GetRadius(); | |
284 | ||
285 | Double_t perp_padrow[2] = {1,0}; //locally in slice | |
286 | ||
287 | Double_t cos_beta = fabs(tangent[0]*perp_padrow[0] + tangent[1]*perp_padrow[1]); | |
288 | return acos(cos_beta); | |
289 | ||
290 | } | |
291 | ||
292 | Bool_t AliL3HoughTrack::GetCrossingPoint(Int_t padrow,Float_t *xyz) | |
293 | { | |
294 | //Assumes the track is given in local coordinates | |
295 | ||
296 | if(!fIsHelix) | |
297 | { | |
298 | printf("AliL3HoughTrack::GetCrossingPoint : Track is not a helix\n"); | |
299 | return 0; | |
300 | } | |
301 | ||
302 | ||
303 | if(!IsLocal()) | |
304 | { | |
305 | printf("GetCrossingPoint: Track is given on global coordinates\n"); | |
306 | return false; | |
307 | } | |
308 | ||
309 | Double_t xHit = fTransform->Row2X(padrow); | |
310 | ||
311 | //xyz[0] = fTransform->Row2X(padrow); | |
312 | xyz[0] = xHit; | |
313 | Double_t aa = (xHit - GetCenterX())*(xHit - GetCenterX()); | |
314 | Double_t r2 = GetRadius()*GetRadius(); | |
315 | if(aa > r2) | |
316 | return false; | |
317 | ||
318 | Double_t aa2 = sqrt(r2 - aa); | |
319 | Double_t y1 = GetCenterY() + aa2; | |
320 | Double_t y2 = GetCenterY() - aa2; | |
321 | xyz[1] = y1; | |
322 | if(fabs(y2) < fabs(y1)) xyz[1] = y2; | |
323 | xyz[2] = 0; //only consider transverse plane | |
324 | ||
325 | return true; | |
326 | } | |
327 | ||
328 | ||
329 | Bool_t AliL3HoughTrack::GetCrossingPoint(Int_t slice,Int_t padrow,Float_t *xyz) | |
330 | { | |
331 | //Calculate the crossing point in transverse plane of this track and given | |
332 | //padrow (y = a*x + b). Point is given in local coordinates in slice. | |
333 | //Assumes the track is given in global coordinates. | |
334 | ||
335 | if(!fIsHelix) | |
336 | { | |
337 | printf("AliL3HoughTrack::GetCrossingPoint : Track is not a helix\n"); | |
338 | return 0; | |
339 | } | |
340 | ||
341 | ||
342 | if(IsLocal()) | |
343 | { | |
344 | printf("GetCrossingPoint: Track is given in local coordintes!!!\n"); | |
345 | return false; | |
346 | } | |
347 | ||
348 | Double_t padrowradii = fTransform->Row2X(padrow); | |
349 | ||
350 | ||
351 | Float_t rotation_angle = (slice*20)*ToRad; | |
352 | ||
353 | Float_t cs,sn; | |
354 | cs = cos(rotation_angle); | |
355 | sn = sin(rotation_angle); | |
356 | ||
357 | Double_t a = -1.*cs/sn; | |
358 | Double_t b = padrowradii/sn; | |
359 | ||
360 | Double_t ycPrime = GetCenterY() - b ; | |
361 | Double_t aa = ( 1. + a * a ) ; | |
362 | Double_t bb = -2. * ( GetCenterX() + a * ycPrime ) ; | |
363 | Double_t cc = ( GetCenterX() * GetCenterX() + ycPrime * ycPrime - GetRadius() * GetRadius() ) ; | |
364 | ||
365 | Double_t racine = bb * bb - 4. * aa * cc ; | |
366 | if ( racine < 0 ) return false ; | |
367 | Double_t rootRacine = sqrt(racine) ; | |
368 | ||
369 | Double_t oneOverA = 1./aa; | |
370 | // | |
371 | // First solution | |
372 | // | |
373 | Double_t x1 = 0.5 * oneOverA * ( -1. * bb + rootRacine ) ; | |
374 | Double_t y1 = a * x1 + b ; | |
375 | Double_t r1 = sqrt(x1*x1+y1*y1); | |
376 | // | |
377 | // Second solution | |
378 | // | |
379 | Double_t x2 = 0.5 * oneOverA * ( -1. * bb - rootRacine ) ; | |
380 | Double_t y2 = a * x2 + b ; | |
381 | Double_t r2 = sqrt(x2*x2+y2*y2); | |
382 | // | |
383 | // Choose close to (0,0) | |
384 | // | |
385 | Double_t xHit ; | |
386 | Double_t yHit ; | |
387 | if ( r1 < r2 ) { | |
388 | xHit = x1 ; | |
389 | yHit = y1 ; | |
390 | } | |
391 | else { | |
392 | xHit = x2 ; | |
393 | yHit = y2 ; | |
394 | } | |
395 | ||
396 | xyz[0] = xHit; | |
397 | xyz[1] = yHit; | |
398 | xyz[2] = 0; | |
399 | ||
400 | fTransform->Global2Local(xyz,slice); | |
401 | ||
402 | return true; | |
403 | } | |
404 | */ |