]>
Commit | Line | Data |
---|---|---|
f8832015 | 1 | #include <TString.h> |
32d38de2 | 2 | #include <TMath.h> |
3 | #include "AliITSUSeed.h" | |
716ccba7 | 4 | #include "AliLog.h" |
3e4e3c23 | 5 | #include "AliESDtrack.h" |
32d38de2 | 6 | using namespace TMath; |
7 | ||
8 | ClassImp(AliITSUSeed) | |
9 | ||
10 | //_________________________________________________________________________ | |
11 | AliITSUSeed::AliITSUSeed() | |
f8832015 | 12 | : fHitsPattern(0) |
b8b59e05 | 13 | ,fNChildren(0) |
f8832015 | 14 | ,fClID(0) |
15 | ,fChi2Glo(0) | |
16 | ,fChi2Cl(0) | |
716ccba7 | 17 | ,fChi2Penalty(0) |
9f42ebf7 | 18 | ,fChi2Match(0) |
19 | ,fChi2ITSSA(0) | |
c61e50c3 | 20 | ,fParent(0) |
08419930 | 21 | #ifdef _ITSU_TUNING_MODE_ // this is used only for tuning histo filling |
22 | ,fOrdBranch(0) | |
23 | ,fOrdCand(0) | |
24 | #endif | |
32d38de2 | 25 | { |
26 | // def c-tor | |
44785f3e | 27 | ResetFMatrix(); |
d9d9b5cd | 28 | for (int i=kNKElem;i--;) fKMatrix[i] = 0.; |
29 | for (int i=kNRElem;i--;) fRMatrix[i] = 0.; | |
30 | fResid[0]=0.; | |
31 | fResid[1]=0.; | |
32 | fCovIYZ[0]=0.; | |
33 | fCovIYZ[1]=0.; | |
34 | fCovIYZ[2]=0.; | |
32d38de2 | 35 | } |
36 | ||
37 | //_________________________________________________________________________ | |
38 | AliITSUSeed::~AliITSUSeed() | |
39 | { | |
716ccba7 | 40 | // d-tor |
32d38de2 | 41 | } |
42 | ||
43 | //_________________________________________________________________________ | |
44 | AliITSUSeed::AliITSUSeed(const AliITSUSeed& src) | |
c61e50c3 | 45 | :AliExternalTrackParam(src) |
f8832015 | 46 | ,fHitsPattern(src.fHitsPattern) |
b8b59e05 | 47 | ,fNChildren(src.fNChildren) |
c61e50c3 | 48 | ,fClID(src.fClID) |
f8832015 | 49 | ,fChi2Glo(src.fChi2Glo) |
50 | ,fChi2Cl(src.fChi2Cl) | |
716ccba7 | 51 | ,fChi2Penalty(src.fChi2Penalty) |
9f42ebf7 | 52 | ,fChi2Match(src.fChi2Match) |
53 | ,fChi2ITSSA(src.fChi2ITSSA) | |
c61e50c3 | 54 | ,fParent(src.fParent) |
08419930 | 55 | #ifdef _ITSU_TUNING_MODE_ // this is used only for tuning histo filling |
56 | ,fOrdBranch(src.fOrdBranch) | |
57 | ,fOrdCand(src.fOrdCand) | |
58 | #endif | |
32d38de2 | 59 | { |
60 | // def c-tor | |
44785f3e | 61 | for (int i=kNFElem;i--;) fFMatrix[i] = src.fFMatrix[i]; |
716ccba7 | 62 | for (int i=kNKElem;i--;) fKMatrix[i] = src.fKMatrix[i]; |
63 | for (int i=kNRElem;i--;) fRMatrix[i] = src.fRMatrix[i]; | |
943e1898 | 64 | fResid[0]=src.fResid[0]; |
65 | fResid[1]=src.fResid[1]; | |
66 | fCovIYZ[0]=src.fCovIYZ[0]; | |
67 | fCovIYZ[1]=src.fCovIYZ[1]; | |
68 | fCovIYZ[2]=src.fCovIYZ[2]; | |
69 | // | |
32d38de2 | 70 | } |
71 | ||
72 | //_________________________________________________________________________ | |
73 | AliITSUSeed &AliITSUSeed::operator=(const AliITSUSeed& src) | |
74 | { | |
75 | // def c-tor | |
76 | if (this == &src) return *this; | |
943e1898 | 77 | this->~AliITSUSeed(); |
78 | new(this) AliITSUSeed(src); | |
32d38de2 | 79 | return *this; |
80 | } | |
f8832015 | 81 | |
82 | //_________________________________________________________________________ | |
83 | void AliITSUSeed::Print(Option_t* opt) const | |
84 | { | |
85 | // print seed info | |
86 | int lr,cl = GetLrCluster(lr); | |
08419930 | 87 | printf("%cLr%d Nchild: %3d Cl:%4d Chi2Glo:%7.2f(%7.2f) Chi2Cl:%7.2f Penalty: %7.2f Mtc:%6.3f Bwd:%6.3f",IsKilled() ? '-':' ', |
9f42ebf7 | 88 | lr,GetNChildren(),cl,GetChi2Glo(),GetChi2GloNrm(),GetChi2Cl(), GetChi2Penalty(), GetChi2ITSTPC(), GetChi2ITSSA()); |
f8832015 | 89 | printf(" |"); |
716ccba7 | 90 | int lrc=0; |
91 | const AliITSUSeed *sdc = this; | |
92 | while(1) { | |
93 | if (lrc<lr) printf("."); | |
94 | else { | |
95 | sdc = sdc->GetParent(lrc); | |
96 | if (!sdc) break; | |
97 | printf("%c",sdc->GetClusterID()<0 ? '.': (sdc->IsFake() ? '-':'+')); | |
98 | } | |
99 | lrc++; | |
100 | } | |
101 | printf("|\n"); | |
f8832015 | 102 | TString opts = opt; opts.ToLower(); |
103 | if (opts.Contains("etp")) AliExternalTrackParam::Print(); | |
104 | if (opts.Contains("parent") && GetParent()) GetParent()->Print(opt); | |
105 | } | |
3dd9c283 | 106 | |
3e4e3c23 | 107 | //______________________________________________________________________________ |
9f42ebf7 | 108 | void AliITSUSeed::InitFromSeed(const AliExternalTrackParam* seed) |
3e4e3c23 | 109 | { |
110 | // init seed from ESD track | |
111 | TObject::Clear(); | |
9f42ebf7 | 112 | AliExternalTrackParam::operator=(*seed); |
3e4e3c23 | 113 | ResetFMatrix(); |
114 | fHitsPattern = 0; | |
115 | fClID = 0; | |
b8b59e05 | 116 | fNChildren = 0; |
3e4e3c23 | 117 | fChi2Glo = fChi2Cl = fChi2Penalty = 0; |
118 | fParent = 0; //!!! | |
119 | } | |
120 | ||
121 | ||
3dd9c283 | 122 | //______________________________________________________________________________ |
123 | Float_t AliITSUSeed::GetChi2GloNrm() const | |
124 | { | |
125 | int ndf = 2*GetNLayersHit() - 5; | |
716ccba7 | 126 | return (ndf>0 ? fChi2Glo/ndf : fChi2Glo) + fChi2Penalty; |
3dd9c283 | 127 | } |
128 | ||
129 | ||
130 | //______________________________________________________________________________ | |
131 | Int_t AliITSUSeed::Compare(const TObject* obj) const | |
132 | { | |
133 | // compare clusters accodring to specific mode | |
f9c7eb32 | 134 | const AliITSUSeed* sd = (const AliITSUSeed*)obj; |
3dd9c283 | 135 | const Float_t kTol = 1e-5; |
136 | if (!IsKilled() && sd->IsKilled()) return -1; | |
137 | if ( IsKilled() &&!sd->IsKilled()) return 1; | |
138 | // | |
716ccba7 | 139 | float chi2This = GetChi2GloNrm(); |
140 | float chi2Other = sd->GetChi2GloNrm(); | |
141 | ||
142 | if (chi2This+kTol<chi2Other) return -1; | |
143 | else if (chi2This-kTol>chi2Other) return 1; | |
3dd9c283 | 144 | return 0; |
145 | } | |
146 | ||
147 | //______________________________________________________________________________ | |
148 | Bool_t AliITSUSeed::IsEqual(const TObject* obj) const | |
149 | { | |
150 | // compare clusters accodring to specific mode | |
f9c7eb32 | 151 | const AliITSUSeed* sd = (const AliITSUSeed*)obj; |
3dd9c283 | 152 | const Float_t kTol = 1e-5; |
153 | if (IsKilled() != sd->IsKilled()) return kFALSE; | |
716ccba7 | 154 | return Abs(GetChi2GloNrm() - sd->GetChi2GloNrm())<kTol; |
3dd9c283 | 155 | } |
44785f3e | 156 | |
157 | //______________________________________________________________________________ | |
158 | Bool_t AliITSUSeed::PropagateToX(Double_t xk, Double_t b) | |
159 | { | |
160 | // Propagate this track to the plane X=xk (cm) in the field "b" (kG) | |
161 | Double_t dx=xk-fX; | |
162 | if (TMath::Abs(dx)<=kAlmost0) return kTRUE; | |
163 | ||
164 | Double_t crv=GetC(b); | |
165 | if (TMath::Abs(b) < kAlmost0Field) crv=0.; | |
166 | Double_t x2r = crv*dx; | |
167 | Double_t f1=fP[2], f2=f1 + x2r; | |
168 | if (TMath::Abs(f1) >= kAlmost1) return kFALSE; | |
169 | if (TMath::Abs(f2) >= kAlmost1) return kFALSE; | |
170 | if (TMath::Abs(fP[4])< kAlmost0) return kFALSE; | |
171 | ||
172 | Double_t &fP0=fP[0], &fP1=fP[1], &fP2=fP[2], &fP3=fP[3], &fP4=fP[4]; | |
173 | Double_t | |
174 | &fC00=fC[0], | |
175 | &fC10=fC[1], &fC11=fC[2], | |
176 | &fC20=fC[3], &fC21=fC[4], &fC22=fC[5], | |
177 | &fC30=fC[6], &fC31=fC[7], &fC32=fC[8], &fC33=fC[9], | |
178 | &fC40=fC[10], &fC41=fC[11], &fC42=fC[12], &fC43=fC[13], &fC44=fC[14]; | |
179 | ||
180 | Double_t r1=TMath::Sqrt((1.-f1)*(1.+f1)), r2=TMath::Sqrt((1.-f2)*(1.+f2)); | |
181 | if (TMath::Abs(r1)<kAlmost0) return kFALSE; | |
182 | if (TMath::Abs(r2)<kAlmost0) return kFALSE; | |
183 | ||
184 | fX=xk; | |
185 | double dy2dx = (f1+f2)/(r1+r2); | |
186 | fP0 += dx*dy2dx; | |
187 | if (TMath::Abs(x2r)<0.05) { | |
188 | fP1 += dx*(r2 + f2*dy2dx)*fP3; // Many thanks to P.Hristov ! | |
189 | fP2 += x2r; | |
190 | } | |
191 | else { | |
192 | // for small dx/R the linear apporximation of the arc by the segment is OK, | |
193 | // but at large dx/R the error is very large and leads to incorrect Z propagation | |
194 | // angle traversed delta = 2*asin(dist_start_end / R / 2), hence the arc is: R*deltaPhi | |
195 | // The dist_start_end is obtained from sqrt(dx^2+dy^2) = x/(r1+r2)*sqrt(2+f1*f2+r1*r2) | |
196 | // Similarly, the rotation angle in linear in dx only for dx<<R | |
197 | double chord = dx*TMath::Sqrt(1+dy2dx*dy2dx); // distance from old position to new one | |
198 | double rot = 2*TMath::ASin(0.5*chord*crv); // angular difference seen from the circle center | |
199 | fP1 += rot/crv*fP3; | |
200 | fP2 = TMath::Sin(rot + TMath::ASin(fP2)); | |
201 | } | |
202 | ||
203 | //f = F - 1 | |
204 | double r1i = 1./r1; | |
205 | double r2i = 1./r2; | |
206 | double tg1 = f1*r1i; | |
207 | double tg2 = f2*r2i; | |
208 | double v0 = 1. + dy2dx*tg2; | |
209 | double v1 = (r1i+r2i)*(dy2dx*(tg1+tg2)+2); | |
210 | double v2 = (r1i+r2i)*v0; | |
211 | // | |
212 | double f24 = dx*crv/fP4; | |
213 | double f02 = dx*v1; | |
214 | double f04 = dx*v2*f24; | |
215 | double f12 = dx*fP3* (f2*v1+dy2dx-tg2); | |
216 | double f13 = dx*r2*v0; | |
217 | double f14 = dx*f24*fP3*(f2*v2+dy2dx-tg2); | |
218 | // | |
219 | //b = C*ft | |
220 | Double_t b00=f02*fC20 + f04*fC40, b01=f12*fC20 + f14*fC40 + f13*fC30; | |
221 | Double_t b02=f24*fC40; | |
222 | Double_t b10=f02*fC21 + f04*fC41, b11=f12*fC21 + f14*fC41 + f13*fC31; | |
223 | Double_t b12=f24*fC41; | |
224 | Double_t b20=f02*fC22 + f04*fC42, b21=f12*fC22 + f14*fC42 + f13*fC32; | |
225 | Double_t b22=f24*fC42; | |
226 | Double_t b40=f02*fC42 + f04*fC44, b41=f12*fC42 + f14*fC44 + f13*fC43; | |
227 | Double_t b42=f24*fC44; | |
228 | Double_t b30=f02*fC32 + f04*fC43, b31=f12*fC32 + f14*fC43 + f13*fC33; | |
229 | Double_t b32=f24*fC43; | |
230 | ||
231 | //a = f*b = f*C*ft | |
232 | Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a02=f02*b22+f04*b42; | |
233 | Double_t a11=f12*b21+f14*b41+f13*b31,a12=f12*b22+f14*b42+f13*b32; | |
234 | Double_t a22=f24*b42; | |
235 | ||
236 | //F*C*Ft = C + (b + bt + a) | |
237 | fC00 += b00 + b00 + a00; | |
238 | fC10 += b10 + b01 + a01; | |
239 | fC20 += b20 + b02 + a02; | |
240 | fC30 += b30; | |
241 | fC40 += b40; | |
242 | fC11 += b11 + b11 + a11; | |
243 | fC21 += b21 + b12 + a12; | |
244 | fC31 += b31; | |
245 | fC41 += b41; | |
246 | fC22 += b22 + b22 + a22; | |
247 | fC32 += b32; | |
248 | fC42 += b42; | |
249 | // | |
250 | // update stored transformation matrix F = Fnew*Fold | |
251 | fFMatrix[kF04] += f04 + f24*fFMatrix[kF02]; | |
252 | fFMatrix[kF14] += f14 + f24*fFMatrix[kF12]; | |
253 | fFMatrix[kF02] += f02; | |
254 | fFMatrix[kF12] += f12; | |
255 | fFMatrix[kF13] += f13; | |
256 | fFMatrix[kF24] += f24; | |
257 | // | |
258 | CheckCovariance(); | |
259 | ||
260 | return kTRUE; | |
261 | } | |
943e1898 | 262 | |
e7d83d38 | 263 | //__________________________________________________________________ |
264 | Int_t AliITSUSeed::GetClusterIndex(Int_t ind) const | |
265 | { | |
266 | // get ind-th cluster index | |
267 | int ncl = 0; | |
268 | const AliITSUSeed* seed = this; | |
269 | while(seed) { | |
c425b166 | 270 | if ( seed->HasCluster() && (ncl++==ind) ) return seed->GetLrClusterID();//GetClusterID(); |
e7d83d38 | 271 | seed = (AliITSUSeed*)seed->GetParent(); |
272 | } | |
273 | return -1; | |
274 | // | |
275 | } | |
276 | ||
716ccba7 | 277 | //______________________________________________________________________________ |
278 | Bool_t AliITSUSeed::RotateToAlpha(Double_t alpha) | |
279 | { | |
280 | // Transform this track to the local coord. system rotated | |
281 | // by angle "alpha" (rad) with respect to the global coord. system. | |
282 | // | |
283 | if (TMath::Abs(fP[2]) >= kAlmost1) { | |
284 | AliError(Form("Precondition is not satisfied: |sin(phi)|>1 ! %f",fP[2])); | |
285 | return kFALSE; | |
286 | } | |
287 | // | |
288 | if (alpha < -TMath::Pi()) alpha += 2*TMath::Pi(); | |
289 | else if (alpha >= TMath::Pi()) alpha -= 2*TMath::Pi(); | |
290 | // | |
291 | Double_t &fP0=fP[0]; | |
292 | Double_t &fP2=fP[2]; | |
293 | Double_t &fC00=fC[0]; | |
294 | Double_t &fC10=fC[1]; | |
295 | Double_t &fC20=fC[3]; | |
296 | Double_t &fC21=fC[4]; | |
297 | Double_t &fC22=fC[5]; | |
298 | Double_t &fC30=fC[6]; | |
299 | Double_t &fC32=fC[8]; | |
300 | Double_t &fC40=fC[10]; | |
301 | Double_t &fC42=fC[12]; | |
302 | // | |
303 | Double_t x=fX; | |
304 | Double_t ca=TMath::Cos(alpha-fAlpha), sa=TMath::Sin(alpha-fAlpha); | |
305 | Double_t sf=fP2, cf=TMath::Sqrt((1.- fP2)*(1.+fP2)); // Improve precision | |
306 | // RS: check if rotation does no invalidate track model (cos(local_phi)>=0, i.e. particle | |
307 | // direction in local frame is along the X axis | |
308 | if ((cf*ca+sf*sa)<0) { | |
309 | AliDebug(1,Form("Rotation failed: local cos(phi) would become %.2f",cf*ca+sf*sa)); | |
310 | return kFALSE; | |
311 | } | |
312 | // | |
313 | Double_t tmp=sf*ca - cf*sa; | |
314 | ||
315 | if (TMath::Abs(tmp) >= kAlmost1) { | |
316 | if (TMath::Abs(tmp) > 1.+ Double_t(FLT_EPSILON)) | |
317 | AliWarning(Form("Rotation failed ! %.10e",tmp)); | |
318 | return kFALSE; | |
319 | } | |
320 | fAlpha = alpha; | |
321 | fX = x*ca + fP0*sa; | |
322 | fP0= -x*sa + fP0*ca; | |
323 | fP2= tmp; | |
324 | ||
325 | if (TMath::Abs(cf)<kAlmost0) { | |
326 | AliError(Form("Too small cosine value %f",cf)); | |
327 | cf = kAlmost0; | |
328 | } | |
329 | ||
330 | Double_t rr=(ca+sf/cf*sa); | |
331 | ||
332 | fC00 *= (ca*ca); | |
333 | fC10 *= ca; | |
334 | fC20 *= ca*rr; | |
335 | fC21 *= rr; | |
336 | fC22 *= rr*rr; | |
337 | fC30 *= ca; | |
338 | fC32 *= rr; | |
339 | fC40 *= ca; | |
340 | fC42 *= rr; | |
341 | // | |
342 | fRMatrix[kR00] = ca; | |
343 | fRMatrix[kR22] = rr; | |
344 | // | |
345 | CheckCovariance(); | |
346 | ||
347 | return kTRUE; | |
348 | } | |
349 | ||
943e1898 | 350 | //______________________________________________________________________________ |
351 | Bool_t AliITSUSeed::GetTrackingXAtXAlpha(double xOther, double alpOther, double bz, double &xdst) | |
352 | { | |
353 | // calculate X and Y in the tracking frame of the track, corresponding to other X,Alpha tracking | |
354 | double ca=TMath::Cos(alpOther-fAlpha), sa=TMath::Sin(alpOther-fAlpha); | |
355 | double &y=fP[0], &sf=fP[2], cf=Sqrt((1.-sf)*(1.+sf)); | |
356 | double eta = xOther - fX*ca - y*sa; | |
357 | double xi = sf*ca - cf*sa; | |
716ccba7 | 358 | if (Abs(xi)>= kAlmost1) return kFALSE; |
943e1898 | 359 | double nu = xi + GetC(bz)*eta; |
716ccba7 | 360 | if (Abs(nu)>= kAlmost1) return kFALSE; |
943e1898 | 361 | xdst = xOther*ca - sa*( y*ca-fX*sa + eta*(xi+nu)/(Sqrt((1.-xi)*(1.+xi)) + Sqrt((1.-nu)*(1.+nu))) ); |
362 | return kTRUE; | |
363 | } | |
364 | ||
365 | //____________________________________________________________________ | |
366 | Double_t AliITSUSeed::GetPredictedChi2(Double_t p[2],Double_t cov[3]) | |
367 | { | |
368 | // Estimate the chi2 of the space point "p" with the cov. matrix "cov" | |
369 | // Store info needed for update and smoothing | |
370 | Double_t sdd = fC[0] + cov[0]; | |
371 | Double_t sdz = fC[1] + cov[1]; | |
372 | Double_t szz = fC[2] + cov[2]; | |
373 | Double_t det = sdd*szz - sdz*sdz; | |
374 | if (TMath::Abs(det) < kAlmost0) return kVeryBig; | |
375 | det = 1./det; | |
376 | fCovIYZ[0] = szz*det; | |
377 | fCovIYZ[1] = -sdz*det; | |
378 | fCovIYZ[2] = sdd*det; | |
379 | double &dy = fResid[0] = p[0] - fP[0]; | |
380 | double &dz = fResid[1] = p[1] - fP[1]; | |
381 | // | |
382 | return dy*(dy*fCovIYZ[0]+dz*fCovIYZ[1]) + dz*(dy*fCovIYZ[1]+dz*(fCovIYZ[2])); | |
383 | // | |
384 | } | |
385 | ||
386 | //____________________________________________________________________ | |
387 | Bool_t AliITSUSeed::Update() | |
388 | { | |
389 | // Update the track parameters with the measurement stored during GetPredictedChi2 | |
390 | // | |
716ccba7 | 391 | Double_t &fP0=fP[0], &fP1=fP[1], &fP2=fP[2], &fP3=fP[3], &fP4=fP[4], |
392 | &fC00=fC[kS00], | |
393 | &fC10=fC[kS10], &fC11=fC[kS11], | |
394 | &fC20=fC[kS20], &fC21=fC[kS21], &fC22=fC[kS22], | |
395 | &fC30=fC[kS30], &fC31=fC[kS31], &fC32=fC[kS32], &fC33=fC[kS33], | |
396 | &fC40=fC[kS40], &fC41=fC[kS41], &fC42=fC[kS42], &fC43=fC[kS43], &fC44=fC[kS44]; | |
943e1898 | 397 | // |
398 | double &r00=fCovIYZ[0],&r01=fCovIYZ[1],&r11=fCovIYZ[2]; | |
399 | double &dy=fResid[0], &dz=fResid[1]; | |
400 | // | |
716ccba7 | 401 | // store info needed for smoothing in the fKMatrix |
402 | double &k00 = fKMatrix[kK00] = fC00*r00+fC10*r01; | |
403 | double &k01 = fKMatrix[kK01] = fC00*r01+fC10*r11; | |
404 | double &k10 = fKMatrix[kK10] = fC10*r00+fC11*r01; | |
405 | double &k11 = fKMatrix[kK11] = fC10*r01+fC11*r11; | |
406 | double &k20 = fKMatrix[kK20] = fC20*r00+fC21*r01; | |
407 | double &k21 = fKMatrix[kK21] = fC20*r01+fC21*r11; | |
408 | double &k30 = fKMatrix[kK30] = fC30*r00+fC31*r01; | |
409 | double &k31 = fKMatrix[kK31] = fC30*r01+fC31*r11; | |
410 | double &k40 = fKMatrix[kK40] = fC40*r00+fC41*r01; | |
411 | double &k41 = fKMatrix[kK41] = fC40*r01+fC41*r11; | |
943e1898 | 412 | // |
413 | Double_t sf=fP2 + k20*dy + k21*dz; | |
414 | if (TMath::Abs(sf) > kAlmost1) return kFALSE; | |
415 | ||
416 | fP0 += k00*dy + k01*dz; | |
417 | fP1 += k10*dy + k11*dz; | |
418 | fP2 = sf; | |
419 | fP3 += k30*dy + k31*dz; | |
420 | fP4 += k40*dy + k41*dz; | |
421 | // | |
422 | Double_t c01=fC10, c02=fC20, c03=fC30, c04=fC40; | |
423 | Double_t c12=fC21, c13=fC31, c14=fC41; | |
424 | ||
425 | fC00-=k00*fC00+k01*fC10; fC10-=k00*c01+k01*fC11; | |
426 | fC20-=k00*c02+k01*c12; fC30-=k00*c03+k01*c13; | |
427 | fC40-=k00*c04+k01*c14; | |
428 | ||
429 | fC11-=k10*c01+k11*fC11; | |
430 | fC21-=k10*c02+k11*c12; fC31-=k10*c03+k11*c13; | |
431 | fC41-=k10*c04+k11*c14; | |
432 | ||
433 | fC22-=k20*c02+k21*c12; fC32-=k20*c03+k21*c13; | |
434 | fC42-=k20*c04+k21*c14; | |
435 | ||
436 | fC33-=k30*c03+k31*c13; | |
437 | fC43-=k30*c04+k31*c14; | |
438 | ||
439 | fC44-=k40*c04+k41*c14; | |
440 | // | |
943e1898 | 441 | CheckCovariance(); |
716ccba7 | 442 | // |
443 | return kTRUE; | |
444 | } | |
445 | ||
446 | ||
447 | //____________________________________________________________________ | |
448 | Bool_t AliITSUSeed::Smooth(Double_t vecL[5],Double_t matL[15]) | |
449 | { | |
450 | // Prepare MBF smoothing auxiliary params for smoothing at prev. point: | |
451 | // \hat{l_N} = 0 | |
452 | // \hat{L_N} = 0 | |
453 | // \tilde{l_j} = -H_j^T N_{j}^{-1} z_j + B_{j}^T \hat{l_j} | |
454 | // \tilde{L_j} = H_j^T N_{j}^{-1} H_j + B_j^T \hat{L_j} B_j | |
455 | // \hat{l_j} = F_j^T \tilde{l_{j+1}} | |
456 | // \hat{L_j} = F_j^T \tilde{L_{j+1}} F_j | |
457 | // | |
458 | // P_{j/N} = P_{j/j} - P_{j/j} \hat{L_j} P_{j/j} | |
459 | // \hat{x_{j/N}} = \hat{x_{j/j}} - P_{j/j} \hat{l_j} | |
460 | // | |
461 | // N^-1 = fCovIYZ | |
462 | // z = fResid | |
463 | // B = I - K H | |
464 | // H = {{1,0,0,0,0},{0,1,0,0,0}} | |
465 | // | |
466 | // calc. \tilde{l_j} | |
467 | // | |
468 | if (GetClusterID()<0) return kTRUE; | |
469 | // | |
470 | ||
471 | double | |
472 | &k00=fKMatrix[kK00],&k01=fKMatrix[kK01], | |
473 | &k10=fKMatrix[kK10],&k11=fKMatrix[kK11], | |
474 | &k20=fKMatrix[kK20],&k21=fKMatrix[kK21], | |
475 | &k30=fKMatrix[kK30],&k31=fKMatrix[kK31], | |
476 | &k40=fKMatrix[kK40],&k41=fKMatrix[kK41]; | |
477 | double | |
478 | &l00=matL[kS00], | |
479 | &l10=matL[kS10], &l11=matL[kS11], | |
480 | &l20=matL[kS20], &l21=matL[kS21], &l22=matL[kS22], | |
481 | &l30=matL[kS30], &l31=matL[kS31], &l32=matL[kS32], &l33=matL[kS33], | |
482 | &l40=matL[kS40], &l41=matL[kS41], &l42=matL[kS42], &l43=matL[kS43], &l44=matL[kS44]; | |
483 | // | |
484 | // calculate correction | |
485 | double corrVec[5]={0},corrMat[15]={0}; | |
486 | corrVec[0] = fC[kS00]*vecL[0] + fC[kS10]*vecL[1] + fC[kS20]*vecL[2] + fC[kS30]*vecL[3] + fC[kS40]*vecL[4]; | |
487 | corrVec[1] = fC[kS10]*vecL[0] + fC[kS11]*vecL[1] + fC[kS21]*vecL[2] + fC[kS31]*vecL[3] + fC[kS41]*vecL[4]; | |
488 | corrVec[2] = fC[kS20]*vecL[0] + fC[kS21]*vecL[1] + fC[kS22]*vecL[2] + fC[kS32]*vecL[3] + fC[kS42]*vecL[4]; | |
489 | corrVec[3] = fC[kS30]*vecL[0] + fC[kS31]*vecL[1] + fC[kS32]*vecL[2] + fC[kS33]*vecL[3] + fC[kS43]*vecL[4]; | |
490 | corrVec[4] = fC[kS40]*vecL[0] + fC[kS41]*vecL[1] + fC[kS42]*vecL[2] + fC[kS43]*vecL[3] + fC[kS44]*vecL[4]; | |
491 | // | |
492 | double *crm = ProdABA(fC,matL); | |
493 | for (int i=0;i<15;i++) corrMat[i] = crm[i]; | |
494 | ||
495 | double vcL0 = vecL[0], vcL1 = vecL[1]; | |
496 | vecL[0] -= k00*vcL0+k10*vcL1+k20*vecL[2]+k30*vecL[3]+k40*vecL[4] + fCovIYZ[0]*fResid[0] + fCovIYZ[1]*fResid[1]; | |
497 | vecL[1] -= k01*vcL0+k11*vcL1+k21*vecL[2]+k31*vecL[3]+k41*vecL[4] + fCovIYZ[1]*fResid[0] + fCovIYZ[2]*fResid[1]; | |
498 | ||
499 | /* | |
500 | double vcL0 = vecL[0], vcL1 = vecL[1]; | |
501 | vecL[0] -= k00*vcL0+k10*vcL1+fKMatrix[kK20]*vecL[2]+k30*vecL[3]+k40*vecL[4] + fCovIYZ[0]*fResid[0] + fCovIYZ[1]*fResid[1]; | |
502 | vecL[1] -= k01*vcL0+fKMatrix[kK11]*vcL1+k21*vecL[2]+k31*vecL[3]+k41*vecL[4] + fCovIYZ[1]*fResid[0] + fCovIYZ[2]*fResid[1]; | |
503 | vecL[3] += fFMatrix[kF13]*vecL[1]; | |
504 | vecL[4] = fFMatrix[kF04]*vecL[0] + fFMatrix[kF14]*vecL[1] + fFMatrix[kF24]*vecL[2] + fFMatrix[kF44]*vecL[4]; | |
505 | vecL[2] += fFMatrix[kF02]*vecL[0] + fFMatrix[kF12]*vecL[1]; | |
506 | // | |
507 | */ | |
508 | // and \hat{l_j} in one go | |
509 | ||
510 | // L = H^T * sg * H + (I-KH)^T * L * (I - KH) | |
511 | double v00 = k00*l00+k10*l10+k20*l20+k30*l30+k40*l40; | |
512 | double v10 = k00*l10+k10*l11+k20*l21+k30*l31+k40*l41; | |
513 | double v20 = k00*l20+k10*l21+k20*l22+k30*l32+k40*l42; | |
514 | double v30 = k00*l30+k10*l31+k20*l32+k30*l33+k40*l43; | |
515 | double v40 = k00*l40+k10*l41+k20*l42+k30*l43+k40*l44; | |
516 | // | |
517 | double v01 = k01*l00+k11*l10+k21*l20+k31*l30+k41*l40; | |
518 | double v11 = k01*l10+k11*l11+k21*l21+k31*l31+k41*l41; | |
519 | double v21 = k01*l20+k11*l21+k21*l22+k31*l32+k41*l42; | |
520 | double v31 = k01*l30+k11*l31+k21*l32+k31*l33+k41*l43; | |
521 | double v41 = k01*l40+k11*l41+k21*l42+k31*l43+k41*l44; | |
522 | // | |
523 | // (H^T * K^T * L * K * H) - (L * K * H) - (H^T * K^T * L) + (H^T*N^-1*H) | |
524 | l00 += k00*v00 + k10*v10 + k20*v20 + k30*v30 + k40*v40 - v00 - v00 + fCovIYZ[0]; | |
525 | l10 += k01*v00 + k11*v10 + k21*v20 + k31*v30 + k41*v40 - v01 - v10 + fCovIYZ[1]; | |
526 | l11 += k01*v01 + k11*v11 + k21*v21 + k31*v31 + k41*v41 - v11 - v11 + fCovIYZ[2]; | |
527 | // | |
528 | l20 -= v20; | |
529 | l21 -= v21; | |
530 | l30 -= v30; | |
531 | l31 -= v31; | |
532 | l40 -= v40; | |
533 | l41 -= v41; | |
534 | // | |
535 | printf("CorrMt:\n"); | |
536 | printf("%+e\n%+e %+e\n%+e %+e %+e\n%+e %+e %+e %+e\n%+e %+e %+e %+e %+e\n", | |
537 | corrMat[kS00],corrMat[kS10],corrMat[kS11],corrMat[kS20],corrMat[kS21],corrMat[kS22], | |
538 | corrMat[kS30],corrMat[kS31],corrMat[kS32],corrMat[kS33], | |
539 | corrMat[kS40],corrMat[kS41],corrMat[kS42],corrMat[kS43],corrMat[kS44]); | |
540 | ||
541 | printf("SMcorr: %+e %+e %+e %+e %+e\n",corrVec[0],corrVec[1],corrVec[2],corrVec[3],corrVec[4]); | |
542 | ||
543 | printf("State : "); this->AliExternalTrackParam::Print(""); | |
544 | // | |
545 | printf("\nBefore transport back (RotElems: %+e %+e)\n",fRMatrix[kR00],fRMatrix[kR22]); | |
546 | printf("Res: %+e %+e | Err: %+e %+e %+e\n",fResid[0],fResid[1],fCovIYZ[0],fCovIYZ[1],fCovIYZ[2]); | |
547 | printf("Lr%d VecL: ",GetLayerID()); for (int i=0;i<5;i++) printf("%+e ",vecL[i]); printf("\n"); | |
548 | // | |
549 | printf("%+e\n%+e %+e\n%+e %+e %+e\n%+e %+e %+e %+e\n%+e %+e %+e %+e %+e\n", | |
550 | matL[kS00],matL[kS10],matL[kS11],matL[kS20],matL[kS21],matL[kS22], | |
551 | matL[kS30],matL[kS31],matL[kS32],matL[kS33],matL[kS40],matL[kS41],matL[kS42],matL[kS43],matL[kS44]); | |
552 | // | |
553 | printf("F: "); for (int i=0;i<kNFElem;i++) printf("%+e ",fFMatrix[i]); printf("\n"); | |
554 | printf("K: "); for (int i=0;i<kNKElem;i++) printf("%+e ",fKMatrix[i]); printf("\n"); | |
555 | // | |
556 | // apply rotation matrix (diagonal) | |
557 | vecL[0] *= fRMatrix[kR00]; | |
558 | vecL[2] *= fRMatrix[kR22]; | |
559 | // | |
560 | l00 *= fRMatrix[kR00]*fRMatrix[kR00]; | |
561 | l10 *= fRMatrix[kR00]; | |
562 | l20 *= fRMatrix[kR22]*fRMatrix[kR00]; | |
563 | l21 *= fRMatrix[kR22]; | |
564 | l22 *= fRMatrix[kR22]*fRMatrix[kR22]; | |
565 | l30 *= fRMatrix[kR00]; | |
566 | l32 *= fRMatrix[kR22]; | |
567 | l40 *= fRMatrix[kR00]; | |
568 | l42 *= fRMatrix[kR22]; | |
569 | // | |
570 | // Apply translation matrix F^T. Note, that fFMatrix keeps non-trivial elems of F-1 = f, except the e-loss coeff f44 | |
571 | // We need F^T*L* F = L + (L*f) + (L*f)^T + f^T * (L*f) | |
572 | // | |
573 | double | |
574 | &f02=fFMatrix[kF02],&f04=fFMatrix[kF04], | |
575 | &f12=fFMatrix[kF12],&f13=fFMatrix[kF13],&f14=fFMatrix[kF14], | |
576 | &f24=fFMatrix[kF24], | |
577 | f44 =fFMatrix[kF44]; | |
578 | // | |
579 | vecL[4] = f04*vecL[0]+f14*vecL[1]+f24*vecL[2]+f44*vecL[4]; | |
580 | vecL[3] += f13*vecL[1]; | |
581 | vecL[2] += f02*vecL[0]+f12*vecL[1]; | |
582 | // | |
583 | f44 -= 1.0; // !!!!! | |
584 | // | |
585 | //b = L*f | |
586 | Double_t b02=l00*f02+l10*f12, b03=l10*f13, b04=l00*f04+l10*f14+l20*f24+l40*f44; | |
587 | Double_t b12=l10*f02+l11*f12, b13=l11*f13, b14=l10*f04+l11*f14+l21*f24+l41*f44; | |
588 | Double_t b22=l20*f02+l21*f12, b23=l21*f13, b24=l20*f04+l21*f14+l22*f24+l42*f44; | |
589 | Double_t b32=l30*f02+l31*f12, b33=l31*f13, b34=l30*f04+l31*f14+l32*f24+l43*f44; | |
590 | Double_t b42=l40*f02+l41*f12, b43=l41*f13, b44=l40*f04+l41*f14+l42*f24+l44*f44; | |
591 | // | |
592 | //a = f^T * b = f^T * L * f, profit from symmetry | |
593 | Double_t a22=f02*b02+f12*b12, a33=f13*b13, a44=f04*b04+f14*b14+f24*b24+f44*b44, | |
594 | a32=f13*b12, //= a23=f02*b03+f12*b13, | |
595 | a42=f02*b04+f12*b14, //f04*b02+f14*b12+f24*b22+f44*b42 = a24 | |
596 | a43=f13*b14; //f04*b03+f14*b13+f24*b23+f44*b43 = a34 | |
597 | // | |
598 | // F^T*L* F = L + (b + b^T + a) | |
599 | l44 += b44 + b44 + a44; | |
600 | l43 += b43 + b34 + a43; | |
601 | l42 += b42 + b24 + a42; | |
602 | l41 += b14; | |
603 | l40 += b04; | |
604 | l33 += b33 + b33 + a33; | |
605 | l32 += b32 + b23 + a32; | |
606 | l31 += b13; | |
607 | l30 += b03; | |
608 | l22 += b22 + b23 + a22; | |
609 | l21 += b12; | |
610 | l20 += b02; | |
611 | // | |
612 | printf("After transport back\n"); | |
613 | printf("Lr%d VecL: ",GetLayerID()); for (int i=0;i<5;i++) printf("%+e ",vecL[i]); printf("\n"); | |
614 | // | |
615 | printf("%+e\n%+e %+e\n%+e %+e %+e\n%+e %+e %+e %+e\n%+e %+e %+e %+e %+e\n", | |
616 | matL[kS00],matL[kS10],matL[kS11],matL[kS20],matL[kS21],matL[kS22], | |
617 | matL[kS30],matL[kS31],matL[kS32],matL[kS33],matL[kS40],matL[kS41],matL[kS42],matL[kS43],matL[kS44]); | |
618 | ||
619 | return kTRUE; | |
620 | } | |
621 | ||
622 | //____________________________________________________________________ | |
623 | Double_t* AliITSUSeed::ProdABA(const double a[15],const double b[15]) const | |
624 | { | |
625 | // product of symmetric matrices A*B*A | |
626 | // | |
627 | const Short_t knd[5][5] = { | |
628 | {kS00,kS10,kS20,kS30,kS40}, | |
629 | {kS10,kS11,kS21,kS31,kS41}, | |
630 | {kS20,kS21,kS22,kS32,kS42}, | |
631 | {kS30,kS31,kS32,kS33,kS43}, | |
632 | {kS40,kS41,kS42,kS43,kS44} | |
633 | }; | |
634 | // | |
635 | static double aba[15]; | |
636 | // 1) ba = B*A | |
637 | double ba[5][5]; | |
638 | for (int i=5;i--;) for (int j=5;j--;) { | |
639 | ba[i][j] = 0; | |
640 | for (int k=5;k--;) ba[i][j] += b[knd[i][k]]*a[knd[k][j]]; | |
641 | } | |
642 | // | |
643 | // 2) A * ba, lower triangle only | |
644 | for (int i=5;i--;) for (int j=i+1;j--;) { | |
645 | aba[knd[i][j]] = 0; | |
646 | for (int k=5;k--;) aba[knd[i][j]] += a[knd[i][k]]*ba[k][j]; | |
647 | } | |
648 | // | |
649 | return &aba[0]; | |
650 | } | |
943e1898 | 651 | |
38997c3c | 652 | //____________________________________________________________________ |
653 | Bool_t AliITSUSeed::ContainsFake() const | |
654 | { | |
655 | // check if the full branch containes a fake cluster | |
656 | const AliITSUSeed* seed = this; | |
657 | while(seed) { | |
658 | if ( seed->IsFake() ) return kTRUE; | |
659 | seed = (AliITSUSeed*)seed->GetParent(); | |
660 | } | |
661 | return kFALSE; | |
662 | } | |
663 | ||
9f42ebf7 | 664 | //__________________________________________________________________ |
665 | Int_t AliITSUSeed::FetchClusterInfo(Int_t *clIDarr) const | |
666 | { | |
667 | // fill cl.id's in the array. The clusters of layer L will be set at slots | |
668 | // clID[2L] (and clID[2L+1] if there is an extra cluster). | |
669 | Int_t lr,ncl=0; | |
670 | const AliITSUSeed* seed = this; | |
671 | do { | |
672 | int clID = seed->GetLrCluster(lr); | |
673 | if (clID>=0) { | |
674 | lr<<=1; | |
675 | clIDarr[ clIDarr[lr]<0 ? lr : lr+1 ] = clID; | |
676 | ncl++; | |
677 | } | |
678 | } while ((seed=(AliITSUSeed*)seed->GetParent())); | |
679 | return ncl; | |
680 | } | |
681 | ||
716ccba7 | 682 | /* |
683 | //____________________________________________________________________ | |
684 | Bool_t AliITSUSeed::Smooth(Double_t vecL[5],Double_t matL[15]) | |
685 | { | |
686 | // Prepare MBF smoothing auxiliary params for smoothing at prev. point: | |
687 | // \hat{l_N} = 0 | |
688 | // \hat{L_N} = 0 | |
689 | // \tilde{l_j} = -H_j^T N_{j}^{-1} z_j + B_{j}^T \hat{l_j} | |
690 | // \tilde{L_j} = H_j^T N_{j}^{-1} H_j + B_j^T \hat{L_j} B_j | |
691 | // \hat{l_j} = F_j^T \tilde{l_{j+1}} | |
692 | // \hat{L_j} = F_j^T \tilde{L_{j+1}} F_j | |
693 | // | |
694 | // P_{j/N} = P_{j/j} - P_{j/j} \hat{L_j} P_{j/j} | |
695 | // \hat{x_{j/N}} = \hat{x_{j/j}} - P_{j/j} \hat{l_j} | |
696 | // | |
697 | // N^-1 = fCovIYZ | |
698 | // z = fResid | |
699 | // B = I - K H | |
700 | // H = {{1,0,0,0,0},{0,1,0,0,0}} | |
701 | // | |
702 | // calc. \tilde{l_j} and \hat{l_j} in one go | |
703 | // | |
704 | if (GetClusterID()<0) return kTRUE; | |
705 | // | |
706 | double | |
707 | &k00=fKMatrix[kK00],&k01=fKMatrix[kK01], | |
708 | &k10=fKMatrix[kK10],&k11=fKMatrix[kK11], | |
709 | &k20=fKMatrix[kK20],&k21=fKMatrix[kK21], | |
710 | &k30=fKMatrix[kK30],&k31=fKMatrix[kK31], | |
711 | &k40=fKMatrix[kK40],&k41=fKMatrix[kK41]; | |
712 | double | |
713 | &matL00=matL[kS00], | |
714 | &matL10=matL[kS01], &matL11=matL[kS11], | |
715 | &matL20=matL[kS20], &matL21=matL[kS21], &matL22=matL[kS22], | |
716 | &matL30=matL[kS30], &matL31=matL[kS31], &matL32=matL[kS32], &matL33=matL[kS33], | |
717 | &matL40=matL[kS40], &matL41=matL[kS41], &matL42=matL[kS42], &matL43=matL[kS43], &matL44=matL[kS44]; | |
718 | // | |
719 | double vcL0 = vecL[0], vcL1 = vecL[1]; | |
720 | vecL[0] -= k00*vcL0+k10*vcL1+fKMatrix[kK20]*vecL[2]+k30*vecL[3]+k40*vecL[4] + fCovIYZ[0]*fResid[0] + fCovIYZ[1]*fResid[1]; | |
721 | vecL[1] -= k01*vcL0+fKMatrix[kK11]*vcL1+k21*vecL[2]+k31*vecL[3]+k41*vecL[4] + fCovIYZ[1]*fResid[0] + fCovIYZ[2]*fResid[1]; | |
722 | vecL[3] += fFMatrix[kF13]*vecL[1]; | |
723 | vecL[4] = fFMatrix[kF04]*vecL[0] + fFMatrix[kF14]*vecL[1] + fFMatrix[kF24]*vecL[2] + fFMatrix[kF44]*vecL[4]; | |
724 | vecL[2] += fFMatrix[kF02]*vecL[0] + fFMatrix[kF12]*vecL[1]; | |
725 | // | |
726 | ||
727 | // L = H^T * sg * H + (I-KH)^T * L * (I - KH) | |
728 | double v00 = k00*matL00+k10*matL10+k20*matL20+k30*matL30+k40*matL40; | |
729 | double v10 = k00*matL10+k10*matL11+k20*matL21+k30*matL31+k40*matL41; | |
730 | double v20 = k00*matL20+k10*matL12+k20*matL22+k30*matL32+k40*matL42; | |
731 | double v30 = k00*matL30+k10*matL13+k20*matL23+k30*matL33+k40*matL43; | |
732 | double v40 = k00*matL40+k10*matL14+k20*matL24+k30*matL34+k40*matL44; | |
733 | // | |
734 | double v01 = k01*matL00+k11*matL10+k21*matL20+k31*matL30+k41*matL40; | |
735 | double v11 = k01*matL01+k11*matL11+k21*matL21+k31*matL31+k41*matL41; | |
736 | double v21 = k01*matL02+k11*matL12+k21*matL22+k31*matL32+k41*matL42; | |
737 | double v31 = k01*matL03+k11*matL13+k21*matL23+k31*matL33+k41*matL43; | |
738 | double v41 = k01*matL04+k11*matL14+k21*matL24+k31*matL34+k41*matL44; | |
739 | // | |
740 | double t00 = k00*matL00+k10*matL01+k20*matL02+k30*matL03+k40*matL04; | |
741 | double t10 = k00*matL10+k10*matL11+k20*matL12+k30*matL13+k40*matL14; | |
742 | double t20 = k00*matL20+k10*matL21+k20*matL22+k30*matL23+k40*matL24; | |
743 | double t30 = k00*matL30+k10*matL31+k20*matL32+k30*matL33+k40*matL34; | |
744 | double t40 = k00*matL40+k10*matL41+k20*matL42+k30*matL43+k40*matL44; | |
745 | // | |
746 | double t01 = k01*matL00+k11*matL01+k21*matL02+k31*matL03+k41*matL04; | |
747 | double t11 = k01*matL10+k11*matL11+k21*matL12+k31*matL13+k41*matL14; | |
748 | double t21 = k01*matL20+k11*matL21+k21*matL22+k31*matL23+k41*matL24; | |
749 | double t31 = k01*matL30+k11*matL31+k21*matL32+k31*matL33+k41*matL34; | |
750 | double t41 = k01*matL40+k11*matL41+k21*matL42+k31*matL43+k41*matL44; | |
751 | // | |
752 | // (H^T * K^T * L * K * H) - (L * K * H) - (H^T * K^T * L) + (H^T*N^-1*H) | |
753 | matL00 += k00*v00+k10*v10+k20*v20*k30*v30+k40*v40 - t00 - v00 + fCovIYZ[0]; | |
754 | matL01 += k01*v00+k11*v10+k21*v20*k31*v30+k41*v40 - t01 - v10 + fCovIYZ[1]; | |
755 | matL10 += k00*v01+k10*v11+k20*v21*k30*v31+k40*v41 - t10 - v01 + fCovIYZ[1]; | |
756 | matL11 += k01*v01+k11*v11+k21*v21*k31*v31+k41*v41 - t11 - v11 + fCovIYZ[2]; | |
757 | // | |
758 | matL20 -= t20; | |
759 | matL21 -= t21; | |
760 | matL30 -= t30; | |
761 | matL31 -= t31; | |
762 | matL40 -= t40; | |
763 | matL41 -= t41; | |
764 | // | |
765 | matL02 -= v20; | |
766 | matL03 -= v30; | |
767 | matL04 -= v40; | |
768 | matL12 -= v21; | |
769 | matL13 -= v31; | |
770 | matL14 -= v41; | |
771 | // | |
772 | printf("Lr%d VecL: ",GetLayerID()); for (int i=0;i<5;i++) printf("%+e ",vecL[i]); printf("\n"); | |
773 | printf("F: "); for (int i=0;i<kNFElem;i++) printf("%+e ",fFMatrix[i]); printf("\n"); | |
774 | printf("K: "); for (int i=0;i<kNKElem;i++) printf("%+e ",fKMatrix[i]); printf("\n"); | |
775 | // | |
776 | for (int j=0;j<5;j++) { | |
777 | for (int i=0;i<5;i++) printf("%+e ",matL[j][i]); printf("\n"); | |
778 | } | |
779 | // | |
943e1898 | 780 | return kTRUE; |
781 | } | |
716ccba7 | 782 | |
783 | */ | |
e7d83d38 | 784 |