Back to home page

Project CMSSW displayed by LXR

 
 

    


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   //  if ( (fabs(invPt)-0.1)/invPtErr > 3.) return 0;
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   std::ostringstream str;
0077     str <<"\t ipt: " << invPt <<"+/-"<<invPtErr
0078         <<"\t pt:  " << pt.value() <<"+/-"<<pt.error()
0079         <<"\t phi: " << phi.value() <<"+/-"<<phi.error()
0080         <<"\t cot: " << cotTheta.value() <<"+/-"<<cotTheta.error()
0081         <<"\t tip: " << tip.value() <<"+/-"<<tip.error()
0082         <<"\t zip: " << zip.value() <<"+/-"<<zip.error()
0083         <<"\t charge: " << charge;
0084   std::cout <<str.str()<< std::endl;
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   // sigma = p0+p1/|pt|;
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   // sigma = p0+p1/|pt|;
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   // sigma = p0+p1/|pt|;
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   // phi = asin(r*rho/2) with asin(x) ~= x+x**3/(2*3)
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   // sigma = p0+p1/pt;
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   // solve equtaion p3*result^3 + p1*result + dphi = 0;
0591   // where:  result = 1/pt
0592   //         dphi = phiL1 - phi0
0593   // Cardan way is used
0594 
0595   double p1, p2, p3;  //parameters p1,p3 are negative by parameter fix in fit!
0596   param(eta, p1, p2, p3);
0597 
0598   double dphi = deltaPhi(phiL1, phi0);  // phi0-phiL1
0599   if (fabs(dphi) < 0.01) {
0600     result = -dphi / p1;
0601   } else {
0602     double q = dphi / 2. / p3;
0603     double p = p1 / 3. / p3;       // positive
0604     double D = q * q + p * p * p;  // positive
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   // pt*sigma(1/pt) = p0+p1*pt;
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 }