File indexing completed on 2024-04-06 12:27:17
0001 #include "RecoMuon/TrackerSeedGenerator/interface/L1MuonPixelTrackFitter.h"
0002
0003 #include "DataFormats/L1GlobalMuonTrigger/interface/L1MuGMTCand.h"
0004
0005 #include "RecoTracker/TkMSParametrization/interface/PixelRecoLineRZ.h"
0006 #include "RecoTracker/TkMSParametrization/interface/PixelRecoUtilities.h"
0007 #include "RecoTracker/PixelTrackFitting/interface/PixelTrackBuilder.h"
0008 #include "DataFormats/GeometryCommonDetAlgo/interface/Measurement1D.h"
0009 #include "MagneticField/Engine/interface/MagneticField.h"
0010 #include "RecoTracker/TkSeedingLayers/interface/SeedingHitSet.h"
0011
0012 #include <cmath>
0013
0014 template <class T>
0015 T sqr(T t) {
0016 return t * t;
0017 }
0018
0019 L1MuonPixelTrackFitter::L1MuonPixelTrackFitter(const edm::ParameterSet& cfg)
0020 : theConfig(cfg),
0021 invPtErrorScale{theConfig.getParameter<double>("invPtErrorScale")},
0022 phiErrorScale{theConfig.getParameter<double>("phiErrorScale")},
0023 cotThetaErrorScale{theConfig.getParameter<double>("cotThetaErrorScale")},
0024 tipErrorScale{theConfig.getParameter<double>("tipErrorScale")},
0025 zipErrorScale{theConfig.getParameter<double>("zipErrorScale")} {}
0026
0027 void L1MuonPixelTrackFitter::setL1Constraint(const L1MuGMTCand& muon) {
0028 thePhiL1 = muon.phiValue() + 0.021817;
0029 theEtaL1 = muon.etaValue();
0030 theChargeL1 = muon.charge();
0031 }
0032
0033 void L1MuonPixelTrackFitter::setPxConstraint(const SeedingHitSet& hits) {
0034 theHit1 = hits[0]->globalPosition();
0035 theHit1 = hits[1]->globalPosition();
0036 }
0037
0038 reco::Track* L1MuonPixelTrackFitter::run(const MagneticField& field,
0039 const std::vector<const TrackingRecHit*>& hits,
0040 const TrackingRegion& region) const {
0041 double phi_vtx_fromHits = (theHit2 - theHit1).phi();
0042
0043 double invPt = valInversePt(phi_vtx_fromHits, thePhiL1, theEtaL1);
0044 double invPtErr = errInversePt(invPt, theEtaL1);
0045
0046 int charge = (invPt > 0.) ? 1 : -1;
0047 double curvature = PixelRecoUtilities::curvature(invPt, field);
0048
0049 Circle circle(theHit1, theHit2, curvature);
0050
0051 double valPhi = this->valPhi(circle, charge);
0052 double valTip = this->valTip(circle, curvature);
0053 double valZip = this->valZip(curvature, theHit1, theHit2);
0054 double valCotTheta = this->valCotTheta(PixelRecoLineRZ(theHit1, theHit2));
0055
0056
0057
0058 PixelTrackBuilder builder;
0059 double valPt = (fabs(invPt) > 1.e-4) ? 1. / fabs(invPt) : 1.e4;
0060 double errPt = invPtErrorScale * invPtErr * valPt * valPt;
0061 double eta = asinh(valCotTheta);
0062 double errPhi = this->errPhi(invPt, eta) * phiErrorScale;
0063 double errCotTheta = this->errCotTheta(invPt, eta) * cotThetaErrorScale;
0064 double errTip = this->errTip(invPt, eta) * tipErrorScale;
0065 double errZip = this->errZip(invPt, eta) * zipErrorScale;
0066
0067 Measurement1D pt(valPt, errPt);
0068 Measurement1D phi(valPhi, errPhi);
0069 Measurement1D cotTheta(valCotTheta, errCotTheta);
0070 Measurement1D tip(valTip, errTip);
0071 Measurement1D zip(valZip, errZip);
0072
0073 float chi2 = 0.;
0074
0075
0076
0077
0078
0079
0080
0081
0082
0083
0084
0085
0086
0087 return builder.build(pt, phi, cotTheta, tip, zip, chi2, charge, hits, &field);
0088 }
0089
0090 double L1MuonPixelTrackFitter::valPhi(const Circle& circle, int charge) const {
0091 Circle::Point zero(0., 0., 0.);
0092 Circle::Vector center = circle.center() - zero;
0093 long double radius = center.perp();
0094 center /= radius;
0095 Circle::Vector dir = center.cross(-charge * Circle::Vector(0., 0., 1.));
0096 return dir.phi();
0097 }
0098
0099 double L1MuonPixelTrackFitter::errPhi(double invPt, double eta) const {
0100
0101
0102
0103 double p0, p1;
0104 int ieta = int(10 * fabs(eta));
0105 switch (ieta) {
0106 case 0: {
0107 p0 = 0.000597506;
0108 p1 = 0.00221057;
0109 break;
0110 }
0111 case 1: {
0112 p0 = 0.000591867;
0113 p1 = 0.00278744;
0114 break;
0115 }
0116 case 2: {
0117 p0 = 0.000635666;
0118 p1 = 0.00207433;
0119 break;
0120 }
0121 case 3: {
0122 p0 = 0.000619086;
0123 p1 = 0.00255121;
0124 break;
0125 }
0126 case 4: {
0127 p0 = 0.000572067;
0128 p1 = 0.00310618;
0129 break;
0130 }
0131 case 5: {
0132 p0 = 0.000596239;
0133 p1 = 0.00288442;
0134 break;
0135 }
0136 case 6: {
0137 p0 = 0.000607608;
0138 p1 = 0.00282996;
0139 break;
0140 }
0141 case 7: {
0142 p0 = 0.000606446;
0143 p1 = 0.00281118;
0144 break;
0145 }
0146 case 8: {
0147 p0 = 0.000594076;
0148 p1 = 0.00280546;
0149 break;
0150 }
0151 case 9: {
0152 p0 = 0.000579615;
0153 p1 = 0.00335534;
0154 break;
0155 }
0156 case 10: {
0157 p0 = 0.000659546;
0158 p1 = 0.00340443;
0159 break;
0160 }
0161 case 11: {
0162 p0 = 0.000659031;
0163 p1 = 0.00343151;
0164 break;
0165 }
0166 case 12: {
0167 p0 = 0.000738391;
0168 p1 = 0.00337297;
0169 break;
0170 }
0171 case 13: {
0172 p0 = 0.000798966;
0173 p1 = 0.00330008;
0174 break;
0175 }
0176 case 14: {
0177 p0 = 0.000702997;
0178 p1 = 0.00562643;
0179 break;
0180 }
0181 case 15: {
0182 p0 = 0.000973417;
0183 p1 = 0.00312666;
0184 break;
0185 }
0186 case 16: {
0187 p0 = 0.000995213;
0188 p1 = 0.00564278;
0189 break;
0190 }
0191 case 17: {
0192 p0 = 0.00121436;
0193 p1 = 0.00572704;
0194 break;
0195 }
0196 case 18: {
0197 p0 = 0.00119216;
0198 p1 = 0.00760204;
0199 break;
0200 }
0201 case 19: {
0202 p0 = 0.00141204;
0203 p1 = 0.0093777;
0204 break;
0205 }
0206 default: {
0207 p0 = 0.00153161;
0208 p1 = 0.00940265;
0209 break;
0210 }
0211 }
0212 return p0 + p1 * fabs(invPt);
0213 }
0214
0215 double L1MuonPixelTrackFitter::valCotTheta(const PixelRecoLineRZ& line) const { return line.cotLine(); }
0216
0217 double L1MuonPixelTrackFitter::errCotTheta(double invPt, double eta) const {
0218
0219
0220
0221 double p0, p1;
0222 int ieta = int(10 * fabs(eta));
0223 switch (ieta) {
0224 case 0: {
0225 p0 = 0.00166115;
0226 p1 = 5.75533e-05;
0227 break;
0228 }
0229 case 1: {
0230 p0 = 0.00157525;
0231 p1 = 0.000707437;
0232 break;
0233 }
0234 case 2: {
0235 p0 = 0.00122246;
0236 p1 = 0.000325456;
0237 break;
0238 }
0239 case 3: {
0240 p0 = 0.000852422;
0241 p1 = 0.000429216;
0242 break;
0243 }
0244 case 4: {
0245 p0 = 0.000637561;
0246 p1 = 0.00122298;
0247 break;
0248 }
0249 case 5: {
0250 p0 = 0.000555766;
0251 p1 = 0.00158096;
0252 break;
0253 }
0254 case 6: {
0255 p0 = 0.000641202;
0256 p1 = 0.00143339;
0257 break;
0258 }
0259 case 7: {
0260 p0 = 0.000803207;
0261 p1 = 0.000648816;
0262 break;
0263 }
0264 case 8: {
0265 p0 = 0.000741394;
0266 p1 = 0.0015289;
0267 break;
0268 }
0269 case 9: {
0270 p0 = 0.000652019;
0271 p1 = 0.00168873;
0272 break;
0273 }
0274 case 10: {
0275 p0 = 0.000716902;
0276 p1 = 0.00257556;
0277 break;
0278 }
0279 case 11: {
0280 p0 = 0.000800409;
0281 p1 = 0.00190563;
0282 break;
0283 }
0284 case 12: {
0285 p0 = 0.000808778;
0286 p1 = 0.00264139;
0287 break;
0288 }
0289 case 13: {
0290 p0 = 0.000775757;
0291 p1 = 0.00318478;
0292 break;
0293 }
0294 case 14: {
0295 p0 = 0.000705781;
0296 p1 = 0.00460576;
0297 break;
0298 }
0299 case 15: {
0300 p0 = 0.000580679;
0301 p1 = 0.00748248;
0302 break;
0303 }
0304 case 16: {
0305 p0 = 0.000561667;
0306 p1 = 0.00767487;
0307 break;
0308 }
0309 case 17: {
0310 p0 = 0.000521626;
0311 p1 = 0.0100178;
0312 break;
0313 }
0314 case 18: {
0315 p0 = 0.00064253;
0316 p1 = 0.0106062;
0317 break;
0318 }
0319 case 19: {
0320 p0 = 0.000636868;
0321 p1 = 0.0140047;
0322 break;
0323 }
0324 default: {
0325 p0 = 0.000682478;
0326 p1 = 0.0163569;
0327 break;
0328 }
0329 }
0330 return p0 + p1 * fabs(invPt);
0331 }
0332
0333 double L1MuonPixelTrackFitter::valTip(const Circle& circle, double curvature) const {
0334 Circle::Point zero(0., 0., 0.);
0335 Circle::Vector center = circle.center() - zero;
0336 long double radius = center.perp();
0337 return radius - 1. / fabs(curvature);
0338 }
0339
0340 double L1MuonPixelTrackFitter::errTip(double invPt, double eta) const {
0341
0342
0343
0344 double p0, p1;
0345 int ieta = int(10 * fabs(eta));
0346 switch (ieta) {
0347 case 0: {
0348 p0 = 0.00392416;
0349 p1 = 0.00551809;
0350 break;
0351 }
0352 case 1: {
0353 p0 = 0.00390391;
0354 p1 = 0.00543244;
0355 break;
0356 }
0357 case 2: {
0358 p0 = 0.0040651;
0359 p1 = 0.00406496;
0360 break;
0361 }
0362 case 3: {
0363 p0 = 0.00387782;
0364 p1 = 0.00797637;
0365 break;
0366 }
0367 case 4: {
0368 p0 = 0.00376798;
0369 p1 = 0.00866894;
0370 break;
0371 }
0372 case 5: {
0373 p0 = 0.0042131;
0374 p1 = 0.00462184;
0375 break;
0376 }
0377 case 6: {
0378 p0 = 0.00392579;
0379 p1 = 0.00784685;
0380 break;
0381 }
0382 case 7: {
0383 p0 = 0.00370472;
0384 p1 = 0.00790174;
0385 break;
0386 }
0387 case 8: {
0388 p0 = 0.00364433;
0389 p1 = 0.00928368;
0390 break;
0391 }
0392 case 9: {
0393 p0 = 0.00387578;
0394 p1 = 0.00640431;
0395 break;
0396 }
0397 case 10: {
0398 p0 = 0.00382464;
0399 p1 = 0.00960763;
0400 break;
0401 }
0402 case 11: {
0403 p0 = 0.0038907;
0404 p1 = 0.0104562;
0405 break;
0406 }
0407 case 12: {
0408 p0 = 0.00392525;
0409 p1 = 0.0106442;
0410 break;
0411 }
0412 case 13: {
0413 p0 = 0.00400634;
0414 p1 = 0.011218;
0415 break;
0416 }
0417 case 14: {
0418 p0 = 0.0036229;
0419 p1 = 0.0156403;
0420 break;
0421 }
0422 case 15: {
0423 p0 = 0.00444317;
0424 p1 = 0.00832987;
0425 break;
0426 }
0427 case 16: {
0428 p0 = 0.00465492;
0429 p1 = 0.0179908;
0430 break;
0431 }
0432 case 17: {
0433 p0 = 0.0049652;
0434 p1 = 0.0216647;
0435 break;
0436 }
0437 case 18: {
0438 p0 = 0.0051395;
0439 p1 = 0.0233692;
0440 break;
0441 }
0442 case 19: {
0443 p0 = 0.0062917;
0444 p1 = 0.0262175;
0445 break;
0446 }
0447 default: {
0448 p0 = 0.00714444;
0449 p1 = 0.0253856;
0450 break;
0451 }
0452 }
0453 return p0 + p1 * fabs(invPt);
0454 }
0455
0456 double L1MuonPixelTrackFitter::valZip(double curv, const GlobalPoint& pinner, const GlobalPoint& pouter) const {
0457
0458
0459
0460 double rho3 = curv * curv * curv;
0461 double r1 = pinner.perp();
0462 double phi1 = r1 * curv / 2 + pinner.perp2() * r1 * rho3 / 48.;
0463 double r2 = pouter.perp();
0464 double phi2 = r2 * curv / 2 + pouter.perp2() * r2 * rho3 / 48.;
0465 double z1 = pinner.z();
0466 double z2 = pouter.z();
0467
0468 return z1 - phi1 / (phi1 - phi2) * (z1 - z2);
0469 }
0470
0471 double L1MuonPixelTrackFitter::errZip(double invPt, double eta) const {
0472
0473
0474
0475 double p0, p1;
0476 int ieta = int(10 * fabs(eta));
0477 switch (ieta) {
0478 case 0: {
0479 p0 = 0.0120743;
0480 p1 = 0;
0481 break;
0482 }
0483 case 1: {
0484 p0 = 0.0110343;
0485 p1 = 0.0051199;
0486 break;
0487 }
0488 case 2: {
0489 p0 = 0.00846487;
0490 p1 = 0.00570084;
0491 break;
0492 }
0493 case 3: {
0494 p0 = 0.00668726;
0495 p1 = 0.00331165;
0496 break;
0497 }
0498 case 4: {
0499 p0 = 0.00467126;
0500 p1 = 0.00578239;
0501 break;
0502 }
0503 case 5: {
0504 p0 = 0.0043042;
0505 p1 = 0.00598517;
0506 break;
0507 }
0508 case 6: {
0509 p0 = 0.00515392;
0510 p1 = 0.00495422;
0511 break;
0512 }
0513 case 7: {
0514 p0 = 0.0060843;
0515 p1 = 0.00320512;
0516 break;
0517 }
0518 case 8: {
0519 p0 = 0.00564942;
0520 p1 = 0.00478876;
0521 break;
0522 }
0523 case 9: {
0524 p0 = 0.00532111;
0525 p1 = 0.0073239;
0526 break;
0527 }
0528 case 10: {
0529 p0 = 0.00579429;
0530 p1 = 0.00952782;
0531 break;
0532 }
0533 case 11: {
0534 p0 = 0.00614229;
0535 p1 = 0.00977795;
0536 break;
0537 }
0538 case 12: {
0539 p0 = 0.00714661;
0540 p1 = 0.00550482;
0541 break;
0542 }
0543 case 13: {
0544 p0 = 0.0066593;
0545 p1 = 0.00999362;
0546 break;
0547 }
0548 case 14: {
0549 p0 = 0.00634922;
0550 p1 = 0.0148156;
0551 break;
0552 }
0553 case 15: {
0554 p0 = 0.00600586;
0555 p1 = 0.0318022;
0556 break;
0557 }
0558 case 16: {
0559 p0 = 0.00676919;
0560 p1 = 0.027456;
0561 break;
0562 }
0563 case 17: {
0564 p0 = 0.00670066;
0565 p1 = 0.0317005;
0566 break;
0567 }
0568 case 18: {
0569 p0 = 0.00752392;
0570 p1 = 0.0347714;
0571 break;
0572 }
0573 case 19: {
0574 p0 = 0.00791425;
0575 p1 = 0.0566665;
0576 break;
0577 }
0578 default: {
0579 p0 = 0.00882372;
0580 p1 = 0.0596858;
0581 break;
0582 }
0583 }
0584 return p0 + p1 * fabs(invPt);
0585 }
0586
0587 double L1MuonPixelTrackFitter::valInversePt(double phi0, double phiL1, double eta) const {
0588 double result = 0.;
0589
0590
0591
0592
0593
0594
0595 double p1, p2, p3;
0596 param(eta, p1, p2, p3);
0597
0598 double dphi = deltaPhi(phiL1, phi0);
0599 if (fabs(dphi) < 0.01) {
0600 result = -dphi / p1;
0601 } else {
0602 double q = dphi / 2. / p3;
0603 double p = p1 / 3. / p3;
0604 double D = q * q + p * p * p;
0605 double u = pow(-q + sqrt(D), 1. / 3.);
0606 double v = -pow(q + sqrt(D), 1. / 3.);
0607 result = u + v;
0608 }
0609 return result;
0610 }
0611
0612 double L1MuonPixelTrackFitter::errInversePt(double invPt, double eta) const {
0613
0614
0615
0616 double p0, p1;
0617 int ieta = int(10 * fabs(eta));
0618 switch (ieta) {
0619 case 0: {
0620 p0 = 0.0196835;
0621 p1 = 0.00517533;
0622 break;
0623 }
0624 case 1: {
0625 p0 = 0.0266583;
0626 p1 = 0.00478101;
0627 break;
0628 }
0629 case 2: {
0630 p0 = 0.0217164;
0631 p1 = 0.00545425;
0632 break;
0633 }
0634 case 3: {
0635 p0 = 0.0197547;
0636 p1 = 0.00552263;
0637 break;
0638 }
0639 case 4: {
0640 p0 = 0.0208778;
0641 p1 = 0.00536009;
0642 break;
0643 }
0644 case 5: {
0645 p0 = 0.024192;
0646 p1 = 0.00521709;
0647 break;
0648 }
0649 case 6: {
0650 p0 = 0.0265315;
0651 p1 = 0.0051897;
0652 break;
0653 }
0654 case 7: {
0655 p0 = 0.0198071;
0656 p1 = 0.00566822;
0657 break;
0658 }
0659 case 8: {
0660 p0 = 0.0361955;
0661 p1 = 0.00486352;
0662 break;
0663 }
0664 case 9: {
0665 p0 = 0.037864;
0666 p1 = 0.00509094;
0667 break;
0668 }
0669 case 10: {
0670 p0 = 0.0382968;
0671 p1 = 0.00612354;
0672 break;
0673 }
0674 case 11: {
0675 p0 = 0.0308326;
0676 p1 = 0.0074234;
0677 break;
0678 }
0679 case 12: {
0680 p0 = 0.0248577;
0681 p1 = 0.00883049;
0682 break;
0683 }
0684 case 13: {
0685 p0 = 0.0279965;
0686 p1 = 0.00888293;
0687 break;
0688 }
0689 case 14: {
0690 p0 = 0.0372582;
0691 p1 = 0.00950252;
0692 break;
0693 }
0694 case 15: {
0695 p0 = 0.0281366;
0696 p1 = 0.0111501;
0697 break;
0698 }
0699 case 16: {
0700 p0 = 0.0421483;
0701 p1 = 0.0109413;
0702 break;
0703 }
0704 case 17: {
0705 p0 = 0.0461798;
0706 p1 = 0.0125824;
0707 break;
0708 }
0709 case 18: {
0710 p0 = 0.0530603;
0711 p1 = 0.0132638;
0712 break;
0713 }
0714 case 19: {
0715 p0 = 0.0601148;
0716 p1 = 0.0147911;
0717 break;
0718 }
0719 default: {
0720 p0 = 0.0552377;
0721 p1 = 0.0155574;
0722 break;
0723 }
0724 }
0725 return p1 + p0 * fabs(invPt);
0726 }
0727
0728 double L1MuonPixelTrackFitter::findPt(double phi0, double phiL1, double eta, int charge) const {
0729 double dphi_min = fabs(deltaPhi(phi0, phiL1));
0730 double pt_best = 1.;
0731 double pt_cur = 1;
0732 while (pt_cur < 10000.) {
0733 double phi_exp = phi0 + getBending(1. / pt_cur, eta, charge);
0734 double dphi = fabs(deltaPhi(phi_exp, phiL1));
0735 if (dphi < dphi_min) {
0736 pt_best = pt_cur;
0737 dphi_min = dphi;
0738 }
0739 if (pt_cur < 10.)
0740 pt_cur += 0.01;
0741 else if (pt_cur < 20.)
0742 pt_cur += 0.025;
0743 else if (pt_cur < 100.)
0744 pt_cur += 0.1;
0745 else
0746 pt_cur += 1;
0747 };
0748 return pt_best;
0749 }
0750
0751 double L1MuonPixelTrackFitter::getBending(double invPt, double eta, int charge) {
0752 double p1, p2, p3;
0753 param(eta, p1, p2, p3);
0754 return charge * p1 * invPt + charge * p2 * invPt * invPt + charge * p3 * invPt * invPt * invPt;
0755 }
0756
0757 double L1MuonPixelTrackFitter::getBendingError(double invPt, double eta) {
0758 int ieta = int(10 * fabs(eta));
0759 double p0, p1;
0760 switch (ieta) {
0761 case 0: {
0762 p0 = 0.0196835;
0763 p1 = 0.00517533;
0764 break;
0765 }
0766 case 1: {
0767 p0 = 0.0266583;
0768 p1 = 0.00478101;
0769 break;
0770 }
0771 case 2: {
0772 p0 = 0.0217164;
0773 p1 = 0.00545425;
0774 break;
0775 }
0776 case 3: {
0777 p0 = 0.0197547;
0778 p1 = 0.00552263;
0779 break;
0780 }
0781 case 4: {
0782 p0 = 0.0208778;
0783 p1 = 0.00536009;
0784 break;
0785 }
0786 case 5: {
0787 p0 = 0.024192;
0788 p1 = 0.00521709;
0789 break;
0790 }
0791 case 6: {
0792 p0 = 0.0265315;
0793 p1 = 0.0051897;
0794 break;
0795 }
0796 case 7: {
0797 p0 = 0.0198071;
0798 p1 = 0.00566822;
0799 break;
0800 }
0801 case 8: {
0802 p0 = 0.0361955;
0803 p1 = 0.00486352;
0804 break;
0805 }
0806 case 9: {
0807 p0 = 0.037864;
0808 p1 = 0.00509094;
0809 break;
0810 }
0811 case 10: {
0812 p0 = 0.0382968;
0813 p1 = 0.00612354;
0814 break;
0815 }
0816 case 11: {
0817 p0 = 0.0308326;
0818 p1 = 0.0074234;
0819 break;
0820 }
0821 case 12: {
0822 p0 = 0.0248577;
0823 p1 = 0.00883049;
0824 break;
0825 }
0826 case 13: {
0827 p0 = 0.0279965;
0828 p1 = 0.00888293;
0829 break;
0830 }
0831 case 14: {
0832 p0 = 0.0372582;
0833 p1 = 0.00950252;
0834 break;
0835 }
0836 case 15: {
0837 p0 = 0.0281366;
0838 p1 = 0.0111501;
0839 break;
0840 }
0841 case 16: {
0842 p0 = 0.0421483;
0843 p1 = 0.0109413;
0844 break;
0845 }
0846 case 17: {
0847 p0 = 0.0461798;
0848 p1 = 0.0125824;
0849 break;
0850 }
0851 case 18: {
0852 p0 = 0.0530603;
0853 p1 = 0.0132638;
0854 break;
0855 }
0856 case 19: {
0857 p0 = 0.0601148;
0858 p1 = 0.0147911;
0859 break;
0860 }
0861 default: {
0862 p0 = 0.0552377;
0863 p1 = 0.0155574;
0864 break;
0865 }
0866 }
0867 return p0 + p1 * sqr(invPt);
0868 }
0869
0870 void L1MuonPixelTrackFitter::param(double eta, double& p1, double& p2, double& p3) {
0871 int ieta = int(10 * fabs(eta));
0872 switch (ieta) {
0873 case 0: {
0874 p1 = -2.68016;
0875 p2 = 0;
0876 p3 = -12.9653;
0877 break;
0878 }
0879 case 1: {
0880 p1 = -2.67864;
0881 p2 = 0;
0882 p3 = -12.0036;
0883 break;
0884 }
0885 case 2: {
0886 p1 = -2.72997;
0887 p2 = 0;
0888 p3 = -10.3468;
0889 break;
0890 }
0891 case 3: {
0892 p1 = -2.68836;
0893 p2 = 0;
0894 p3 = -12.3369;
0895 break;
0896 }
0897 case 4: {
0898 p1 = -2.66885;
0899 p2 = 0;
0900 p3 = -11.589;
0901 break;
0902 }
0903 case 5: {
0904 p1 = -2.64932;
0905 p2 = 0;
0906 p3 = -12.7176;
0907 break;
0908 }
0909 case 6: {
0910 p1 = -2.64513;
0911 p2 = 0;
0912 p3 = -11.3912;
0913 break;
0914 }
0915 case 7: {
0916 p1 = -2.64487;
0917 p2 = 0;
0918 p3 = -9.83239;
0919 break;
0920 }
0921 case 8: {
0922 p1 = -2.58875;
0923 p2 = 0;
0924 p3 = -10.0541;
0925 break;
0926 }
0927 case 9: {
0928 p1 = -2.5551;
0929 p2 = 0;
0930 p3 = -9.75757;
0931 break;
0932 }
0933 case 10: {
0934 p1 = -2.55294;
0935 p2 = 0;
0936 p3 = -7.25238;
0937 break;
0938 }
0939 case 11: {
0940 p1 = -2.45333;
0941 p2 = 0;
0942 p3 = -7.966;
0943 break;
0944 }
0945 case 12: {
0946 p1 = -2.29283;
0947 p2 = 0;
0948 p3 = -7.03231;
0949 break;
0950 }
0951 case 13: {
0952 p1 = -2.13167;
0953 p2 = 0;
0954 p3 = -4.29182;
0955 break;
0956 }
0957 case 14: {
0958 p1 = -1.9102;
0959 p2 = 0;
0960 p3 = -3.21295;
0961 break;
0962 }
0963 case 15: {
0964 p1 = -1.74552;
0965 p2 = 0;
0966 p3 = -0.531827;
0967 break;
0968 }
0969 case 16: {
0970 p1 = -1.53642;
0971 p2 = 0;
0972 p3 = -1.74057;
0973 break;
0974 }
0975 case 17: {
0976 p1 = -1.39446;
0977 p2 = 0;
0978 p3 = -1.56819;
0979 break;
0980 }
0981 case 18: {
0982 p1 = -1.26176;
0983 p2 = 0;
0984 p3 = -0.598631;
0985 break;
0986 }
0987 case 19: {
0988 p1 = -1.14133;
0989 p2 = 0;
0990 p3 = -0.0941055;
0991 break;
0992 }
0993 default: {
0994 p1 = -1.02086;
0995 p2 = 0;
0996 p3 = -0.100491;
0997 break;
0998 }
0999 }
1000 }
1001
1002 double L1MuonPixelTrackFitter::deltaPhi(double phi1, double phi2) const {
1003 while (phi1 >= 2 * M_PI)
1004 phi1 -= 2 * M_PI;
1005 while (phi2 >= 2 * M_PI)
1006 phi2 -= 2 * M_PI;
1007 while (phi1 < 0)
1008 phi1 += 2 * M_PI;
1009 while (phi2 < 0)
1010 phi2 += 2 * M_PI;
1011 double dPhi = phi2 - phi1;
1012
1013 if (dPhi > M_PI)
1014 dPhi -= 2 * M_PI;
1015 if (dPhi < -M_PI)
1016 dPhi += 2 * M_PI;
1017
1018 return dPhi;
1019 }