3 #include "AliITSUSeed.h"
8 //_________________________________________________________________________
9 AliITSUSeed::AliITSUSeed()
20 //_________________________________________________________________________
21 AliITSUSeed::~AliITSUSeed()
26 //_________________________________________________________________________
27 AliITSUSeed::AliITSUSeed(const AliITSUSeed& src)
28 :AliExternalTrackParam(src)
29 ,fHitsPattern(src.fHitsPattern)
31 ,fChi2Glo(src.fChi2Glo)
36 for (int i=kNFElem;i--;) fFMatrix[i] = src.fFMatrix[i];
37 for (int i=kNBElem;i--;) fBMatrix[i] = src.fBMatrix[i];
38 fResid[0]=src.fResid[0];
39 fResid[1]=src.fResid[1];
40 fCovIYZ[0]=src.fCovIYZ[0];
41 fCovIYZ[1]=src.fCovIYZ[1];
42 fCovIYZ[2]=src.fCovIYZ[2];
46 //_________________________________________________________________________
47 AliITSUSeed &AliITSUSeed::operator=(const AliITSUSeed& src)
50 if (this == &src) return *this;
52 new(this) AliITSUSeed(src);
56 //_________________________________________________________________________
57 void AliITSUSeed::Print(Option_t* opt) const
60 int lr,cl = GetLrCluster(lr);
61 printf("%cLr%d Cl:%4d Chi2Glo:%7.2f(%7.2f) Chi2Cl:",IsKilled() ? '-':' ',
62 lr,cl,GetChi2Glo(),GetChi2GloNrm());
63 cl<0 ? printf(" NA ") : printf("%7.2f",GetChi2Cl());
65 for (int i=0;i<=12;i++) printf("%c",HasClusterOnLayer(i) ? '+':'-'); printf("|\n");
66 TString opts = opt; opts.ToLower();
67 if (opts.Contains("etp")) AliExternalTrackParam::Print();
68 if (opts.Contains("parent") && GetParent()) GetParent()->Print(opt);
71 //______________________________________________________________________________
72 Float_t AliITSUSeed::GetChi2GloNrm() const
74 int ndf = 2*GetNLayersHit() - 5;
75 return ndf>0 ? fChi2Glo/ndf : fChi2Glo;
79 //______________________________________________________________________________
80 Int_t AliITSUSeed::Compare(const TObject* obj) const
82 // compare clusters accodring to specific mode
83 const AliITSUSeed* sd = (const AliITSUSeed*)obj;
84 const Float_t kTol = 1e-5;
85 if (!IsKilled() && sd->IsKilled()) return -1;
86 if ( IsKilled() &&!sd->IsKilled()) return 1;
88 if (GetChi2Glo()+kTol<sd->GetChi2Glo()) return -1;
89 else if (GetChi2Glo()-kTol>sd->GetChi2Glo()) return 1;
93 //______________________________________________________________________________
94 Bool_t AliITSUSeed::IsEqual(const TObject* obj) const
96 // compare clusters accodring to specific mode
97 const AliITSUSeed* sd = (const AliITSUSeed*)obj;
98 const Float_t kTol = 1e-5;
99 if (IsKilled() != sd->IsKilled()) return kFALSE;
100 return Abs(GetChi2Glo() - sd->GetChi2Glo())<kTol;
103 //______________________________________________________________________________
104 Bool_t AliITSUSeed::PropagateToX(Double_t xk, Double_t b)
106 // Propagate this track to the plane X=xk (cm) in the field "b" (kG)
108 if (TMath::Abs(dx)<=kAlmost0) return kTRUE;
110 Double_t crv=GetC(b);
111 if (TMath::Abs(b) < kAlmost0Field) crv=0.;
112 Double_t x2r = crv*dx;
113 Double_t f1=fP[2], f2=f1 + x2r;
114 if (TMath::Abs(f1) >= kAlmost1) return kFALSE;
115 if (TMath::Abs(f2) >= kAlmost1) return kFALSE;
116 if (TMath::Abs(fP[4])< kAlmost0) return kFALSE;
118 Double_t &fP0=fP[0], &fP1=fP[1], &fP2=fP[2], &fP3=fP[3], &fP4=fP[4];
121 &fC10=fC[1], &fC11=fC[2],
122 &fC20=fC[3], &fC21=fC[4], &fC22=fC[5],
123 &fC30=fC[6], &fC31=fC[7], &fC32=fC[8], &fC33=fC[9],
124 &fC40=fC[10], &fC41=fC[11], &fC42=fC[12], &fC43=fC[13], &fC44=fC[14];
126 Double_t r1=TMath::Sqrt((1.-f1)*(1.+f1)), r2=TMath::Sqrt((1.-f2)*(1.+f2));
127 if (TMath::Abs(r1)<kAlmost0) return kFALSE;
128 if (TMath::Abs(r2)<kAlmost0) return kFALSE;
131 double dy2dx = (f1+f2)/(r1+r2);
133 if (TMath::Abs(x2r)<0.05) {
134 fP1 += dx*(r2 + f2*dy2dx)*fP3; // Many thanks to P.Hristov !
138 // for small dx/R the linear apporximation of the arc by the segment is OK,
139 // but at large dx/R the error is very large and leads to incorrect Z propagation
140 // angle traversed delta = 2*asin(dist_start_end / R / 2), hence the arc is: R*deltaPhi
141 // The dist_start_end is obtained from sqrt(dx^2+dy^2) = x/(r1+r2)*sqrt(2+f1*f2+r1*r2)
142 // Similarly, the rotation angle in linear in dx only for dx<<R
143 double chord = dx*TMath::Sqrt(1+dy2dx*dy2dx); // distance from old position to new one
144 double rot = 2*TMath::ASin(0.5*chord*crv); // angular difference seen from the circle center
146 fP2 = TMath::Sin(rot + TMath::ASin(fP2));
154 double v0 = 1. + dy2dx*tg2;
155 double v1 = (r1i+r2i)*(dy2dx*(tg1+tg2)+2);
156 double v2 = (r1i+r2i)*v0;
158 double f24 = dx*crv/fP4;
160 double f04 = dx*v2*f24;
161 double f12 = dx*fP3* (f2*v1+dy2dx-tg2);
162 double f13 = dx*r2*v0;
163 double f14 = dx*f24*fP3*(f2*v2+dy2dx-tg2);
166 Double_t b00=f02*fC20 + f04*fC40, b01=f12*fC20 + f14*fC40 + f13*fC30;
167 Double_t b02=f24*fC40;
168 Double_t b10=f02*fC21 + f04*fC41, b11=f12*fC21 + f14*fC41 + f13*fC31;
169 Double_t b12=f24*fC41;
170 Double_t b20=f02*fC22 + f04*fC42, b21=f12*fC22 + f14*fC42 + f13*fC32;
171 Double_t b22=f24*fC42;
172 Double_t b40=f02*fC42 + f04*fC44, b41=f12*fC42 + f14*fC44 + f13*fC43;
173 Double_t b42=f24*fC44;
174 Double_t b30=f02*fC32 + f04*fC43, b31=f12*fC32 + f14*fC43 + f13*fC33;
175 Double_t b32=f24*fC43;
178 Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a02=f02*b22+f04*b42;
179 Double_t a11=f12*b21+f14*b41+f13*b31,a12=f12*b22+f14*b42+f13*b32;
180 Double_t a22=f24*b42;
182 //F*C*Ft = C + (b + bt + a)
183 fC00 += b00 + b00 + a00;
184 fC10 += b10 + b01 + a01;
185 fC20 += b20 + b02 + a02;
188 fC11 += b11 + b11 + a11;
189 fC21 += b21 + b12 + a12;
192 fC22 += b22 + b22 + a22;
196 // update stored transformation matrix F = Fnew*Fold
197 fFMatrix[kF04] += f04 + f24*fFMatrix[kF02];
198 fFMatrix[kF14] += f14 + f24*fFMatrix[kF12];
199 fFMatrix[kF02] += f02;
200 fFMatrix[kF12] += f12;
201 fFMatrix[kF13] += f13;
202 fFMatrix[kF24] += f24;
209 //______________________________________________________________________________
210 Bool_t AliITSUSeed::GetTrackingXAtXAlpha(double xOther, double alpOther, double bz, double &xdst)
212 // calculate X and Y in the tracking frame of the track, corresponding to other X,Alpha tracking
213 double ca=TMath::Cos(alpOther-fAlpha), sa=TMath::Sin(alpOther-fAlpha);
214 double &y=fP[0], &sf=fP[2], cf=Sqrt((1.-sf)*(1.+sf));
215 double eta = xOther - fX*ca - y*sa;
216 double xi = sf*ca - cf*sa;
217 if (xi>= kAlmost1) return kFALSE;
218 double nu = xi + GetC(bz)*eta;
219 if (nu>= kAlmost1) return kFALSE;
220 xdst = xOther*ca - sa*( y*ca-fX*sa + eta*(xi+nu)/(Sqrt((1.-xi)*(1.+xi)) + Sqrt((1.-nu)*(1.+nu))) );
224 //____________________________________________________________________
225 Double_t AliITSUSeed::GetPredictedChi2(Double_t p[2],Double_t cov[3])
227 // Estimate the chi2 of the space point "p" with the cov. matrix "cov"
228 // Store info needed for update and smoothing
229 Double_t sdd = fC[0] + cov[0];
230 Double_t sdz = fC[1] + cov[1];
231 Double_t szz = fC[2] + cov[2];
232 Double_t det = sdd*szz - sdz*sdz;
233 if (TMath::Abs(det) < kAlmost0) return kVeryBig;
235 fCovIYZ[0] = szz*det;
236 fCovIYZ[1] = -sdz*det;
237 fCovIYZ[2] = sdd*det;
238 double &dy = fResid[0] = p[0] - fP[0];
239 double &dz = fResid[1] = p[1] - fP[1];
241 return dy*(dy*fCovIYZ[0]+dz*fCovIYZ[1]) + dz*(dy*fCovIYZ[1]+dz*(fCovIYZ[2]));
245 //____________________________________________________________________
246 Bool_t AliITSUSeed::Update()
248 // Update the track parameters with the measurement stored during GetPredictedChi2
250 Double_t &fP0=fP[0], &fP1=fP[1], &fP2=fP[2], &fP3=fP[3], &fP4=fP[4];
253 &fC10=fC[1], &fC11=fC[2],
254 &fC20=fC[3], &fC21=fC[4], &fC22=fC[5],
255 &fC30=fC[6], &fC31=fC[7], &fC32=fC[8], &fC33=fC[9],
256 &fC40=fC[10], &fC41=fC[11], &fC42=fC[12], &fC43=fC[13], &fC44=fC[14];
258 double &r00=fCovIYZ[0],&r01=fCovIYZ[1],&r11=fCovIYZ[2];
259 double &dy=fResid[0], &dz=fResid[1];
261 // store info needed for smoothing in the fBMatrix
262 double &k00 = fBMatrix[kB00] = fC00*r00+fC10*r01;
263 double &k01 = fBMatrix[kB01] = fC00*r01+fC10*r11;
264 double &k10 = fBMatrix[kB10] = fC10*r00+fC11*r01;
265 double &k11 = fBMatrix[kB11] = fC10*r01+fC11*r11;
266 double &k20 = fBMatrix[kB20] = fC20*r00+fC21*r01;
267 double &k21 = fBMatrix[kB21] = fC20*r01+fC21*r11;
268 double &k30 = fBMatrix[kB30] = fC30*r00+fC31*r01;
269 double &k31 = fBMatrix[kB31] = fC30*r01+fC31*r11;
270 double &k40 = fBMatrix[kB40] = fC40*r00+fC41*r01;
271 double &k41 = fBMatrix[kB41] = fC40*r01+fC41*r11;
273 Double_t sf=fP2 + k20*dy + k21*dz;
274 if (TMath::Abs(sf) > kAlmost1) return kFALSE;
276 fP0 += k00*dy + k01*dz;
277 fP1 += k10*dy + k11*dz;
279 fP3 += k30*dy + k31*dz;
280 fP4 += k40*dy + k41*dz;
282 Double_t c01=fC10, c02=fC20, c03=fC30, c04=fC40;
283 Double_t c12=fC21, c13=fC31, c14=fC41;
285 fC00-=k00*fC00+k01*fC10; fC10-=k00*c01+k01*fC11;
286 fC20-=k00*c02+k01*c12; fC30-=k00*c03+k01*c13;
287 fC40-=k00*c04+k01*c14;
289 fC11-=k10*c01+k11*fC11;
290 fC21-=k10*c02+k11*c12; fC31-=k10*c03+k11*c13;
291 fC41-=k10*c04+k11*c14;
293 fC22-=k20*c02+k21*c12; fC32-=k20*c03+k21*c13;
294 fC42-=k20*c04+k21*c14;
296 fC33-=k30*c03+k31*c13;
297 fC43-=k30*c04+k31*c14;
299 fC44-=k40*c04+k41*c14;