63 if (listHitSector.size() < 2) {
70 HepSymMatrix error(2, 0);
71 HepVector gloEta(2, 0);
72 HepSymMatrix gloError(2, 0);
76 HepVector globalPar(4, 0);
77 HepSymMatrix globalErr(4, 0);
80 globalPar.sub(1, eta);
81 globalErr.sub(1, error);
85 globalPar.sub(3, eta);
86 globalErr.sub(3, error);
120 double x_mea = hit->getPositionX();
121 double y_mea = hit->getPositionY();
122 double z_mea = hit->getPositionZ();
125 double x_pre, y_pre, z_pre, dx, dy, dz;
129 HepMatrix errors(2, 2, 0);
130 HepMatrix A(2, 4, 0);
134 double hit_xErr;
double hit_yErr;
double hit_zErr;
135 HepMatrix globalHitErr(3, 3, 0);
139 A[ MY ][ BY ] = x_mea;
141 A[ MZ ][ BZ ] = x_mea;
162 int nStrips = hit->getPhiStripMax() - hit->getPhiStripMin() + 1;
163 double dn = nStrips - 1.5;
164 double factor = std::pow((0.9 + 0.4 * dn * dn), 1.5) * 0.60;
165 hit_localPhiErr = hit_localPhiErr *
sqrt(factor);
167 nStrips = hit->getZStripMax() - hit->getZStripMin() + 1;
169 factor = std::pow((0.9 + 0.4 * dn * dn), 1.5) * 0.55;
170 hit_localZErr = hit_localZErr *
sqrt(factor);
174 double sinphi = globalOrigin[1] / globalOrigin.mag();
175 double cosphi = globalOrigin[0] / globalOrigin.mag();
177 hit_xErr = hit_localPhiErr * sinphi;
178 hit_yErr = hit_localPhiErr * cosphi;
179 hit_zErr = hit_localZErr;
182 globalHitErr[0][0] = pow(hit_xErr, 2);
183 globalHitErr[0][1] = (hit_xErr) * (hit_yErr);
184 globalHitErr[0][2] = 0;
185 globalHitErr[1][1] = pow(hit_yErr, 2);;
186 globalHitErr[1][0] = (hit_xErr) * (hit_yErr);
187 globalHitErr[1][2] = 0;
188 globalHitErr[2][2] = pow(hit_zErr, 2);
189 globalHitErr[2][0] = 0;
190 globalHitErr[2][1] = 0;
207 (hit->getXStripMax() - hit->getXStripMin() + 1) /
sqrt(12);
209 (hit->getYStripMax() - hit->getYStripMin() + 1) /
sqrt(12);
212 globalHitErr[0][0] = pow(hit_xErr, 2);
213 globalHitErr[0][1] = 0.;
214 globalHitErr[0][2] = 0.;
215 globalHitErr[1][1] = pow(hit_yErr, 2);
216 globalHitErr[1][0] = 0.;
217 globalHitErr[1][2] = 0.;
218 globalHitErr[2][2] = pow(hit_zErr, 2);;
219 globalHitErr[2][0] = 0.;
220 globalHitErr[2][1] = 0.;
226 B2FATAL(
"In KLMTrackFitter globalDistanceToHit, hit is neither from E/B-KLM.");
231 double distance =
sqrt(dx * dx + dy * dy + dz * dz);
234 error =
sqrt(errors[ MY ][ MY ] +
243 sigma = distance / error;
256 int depDir,
int indDir)
261 HepMatrix globalHitErr(3, 3, 0);
267 int noPoints = hitList.size();
270 HepMatrix A(noPoints, 2, 0);
273 HepVector y(noPoints, 0);
277 HepDiagMatrix V_y_inverse(noPoints, 0);
280 HepSymMatrix V_A, V_A_inverse;
284 for (std::list< KLMHit2d* >::iterator iHit = hitList.begin(); iHit != hitList.end(); ++iHit) {
288 CLHEP::Hep3Vector globalPos(hit->getPositionX(), hit->getPositionY(), hit->getPositionZ());
300 int nStrips = hit->getPhiStripMax() - hit->getPhiStripMin() + 1;
301 double dn = nStrips - 1.5;
302 double factor = std::pow((0.9 + 0.4 * dn * dn), 1.5) * 0.60;
303 hit_localPhiErr = hit_localPhiErr *
sqrt(factor);
305 nStrips = hit->getZStripMax() - hit->getZStripMin() + 1;
307 factor = std::pow((0.9 + 0.4 * dn * dn), 1.5) * 0.55;
308 hit_localZErr = hit_localZErr *
sqrt(factor);
312 double sinphi = globalOrigin[1] / globalOrigin.mag();
313 double cosphi = globalOrigin[0] / globalOrigin.mag();
315 globalHitErr[0][0] = pow(hit_localPhiErr * sinphi, 2);
316 globalHitErr[0][1] = (hit_localPhiErr * sinphi) * (hit_localPhiErr * cosphi);
317 globalHitErr[0][2] = 0.;
318 globalHitErr[1][1] = pow(hit_localPhiErr * cosphi, 2);
319 globalHitErr[1][0] = (hit_localPhiErr * sinphi) * (hit_localPhiErr * cosphi);
320 globalHitErr[1][2] = 0.;
321 globalHitErr[2][2] = pow(hit_localZErr, 2);
322 globalHitErr[2][0] = 0.;
323 globalHitErr[2][1] = 0.;
328 indPos = globalPos.x();
332 indPos = globalPos.y();
336 indPos = globalPos.z();
340 B2DEBUG(20,
"error in bklm in-TrackFit: illegal direction");
347 depPos = globalPos.x();
351 depPos = globalPos.y();
355 depPos = globalPos.z();
359 B2DEBUG(20,
"error in bklm dep-TrackFit: illegal direction");
369 (hit->getXStripMax() - hit->getXStripMin() + 1) /
sqrt(12);
371 (hit->getYStripMax() - hit->getYStripMin() + 1) /
sqrt(12);
375 12) <<
" Vec_y: " << hit_yErr *
sqrt(12));
377 globalHitErr[0][0] = pow(hit_xErr, 2);
378 globalHitErr[0][1] = 0.;
379 globalHitErr[0][2] = 0.;
380 globalHitErr[1][1] = pow(hit_yErr, 2);
381 globalHitErr[1][0] = 0.;
382 globalHitErr[1][2] = 0.;
383 globalHitErr[2][2] = pow(hit_zErr, 2);
384 globalHitErr[2][0] = 0.;
385 globalHitErr[2][1] = 0.;
392 indPos = globalPos.x();
396 indPos = globalPos.y();
400 indPos = globalPos.z();
404 B2DEBUG(20,
"error in EKLM in-TrackFit: illegal direction");
411 depPos = globalPos.x();
415 depPos = globalPos.y();
419 depPos = globalPos.z();
423 B2DEBUG(20,
"error in EKLM dep-TrackFit: illegal direction");
429 A[ n ][ 1 ] = indPos;
433 double error_raw = globalHitErr[indDir][indDir] + globalHitErr[depDir][depDir];
434 double weight = error_raw;
436 V_y_inverse[ n ][ n ] = 1.0 / weight;
438 V_y_inverse[ n ][ n ] = DBL_MAX;
444 V_A_inverse = V_y_inverse.similarityT(A);
447 V_A = V_A_inverse.inverse(ierr);
449 eta = V_A * A.T() * V_y_inverse * y;
453 HepMatrix residual = y - A * eta;
454 HepMatrix chisqr = residual.T() * V_y_inverse * residual;
456 return (chisqr.trace());