Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:22:03

0001 #include "L1Trigger/TrackFindingTracklet/interface/TrackletCalculatorDisplaced.h"
0002 #include "L1Trigger/TrackFindingTracklet/interface/Settings.h"
0003 #include "L1Trigger/TrackFindingTracklet/interface/Globals.h"
0004 #include "L1Trigger/TrackFindingTracklet/interface/Tracklet.h"
0005 #include "L1Trigger/TrackFindingTracklet/interface/Stub.h"
0006 #include "L1Trigger/TrackFindingTracklet/interface/L1TStub.h"
0007 #include "L1Trigger/TrackFindingTracklet/interface/IMATH_TrackletCalculator.h"
0008 
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010 #include "FWCore/Utilities/interface/Exception.h"
0011 #include "DataFormats/Math/interface/deltaPhi.h"
0012 
0013 using namespace std;
0014 using namespace trklet;
0015 
0016 TrackletCalculatorDisplaced::TrackletCalculatorDisplaced(string name, Settings const& settings, Globals* global)
0017     : ProcessBase(name, settings, global) {
0018   for (unsigned int ilayer = 0; ilayer < N_LAYER; ilayer++) {
0019     vector<TrackletProjectionsMemory*> tmp(settings.nallstubs(ilayer), nullptr);
0020     trackletprojlayers_.push_back(tmp);
0021   }
0022 
0023   for (unsigned int idisk = 0; idisk < N_DISK; idisk++) {
0024     vector<TrackletProjectionsMemory*> tmp(settings.nallstubs(idisk + N_LAYER), nullptr);
0025     trackletprojdisks_.push_back(tmp);
0026   }
0027 
0028   layer_ = 0;
0029   disk_ = 0;
0030 
0031   string name1 = name.substr(1);  //this is to correct for "TCD" having one more letter then "TC"
0032   if (name1[3] == 'L')
0033     layer_ = name1[4] - '0';
0034   if (name1[3] == 'D')
0035     disk_ = name1[4] - '0';
0036 
0037   // set TC index
0038   iSeed_ = 0;
0039 
0040   int iTC = name1[9] - 'A';
0041 
0042   if (name1.substr(3, 6) == "L3L4L2")
0043     iSeed_ = 8;
0044   else if (name1.substr(3, 6) == "L5L6L4")
0045     iSeed_ = 9;
0046   else if (name1.substr(3, 6) == "L2L3D1")
0047     iSeed_ = 10;
0048   else if (name1.substr(3, 6) == "D1D2L2")
0049     iSeed_ = 11;
0050 
0051   assert(iSeed_ != 0);
0052 
0053   TCIndex_ = (iSeed_ << 4) + iTC;
0054   assert(TCIndex_ >= 128 && TCIndex_ < 191);
0055 
0056   assert((layer_ != 0) || (disk_ != 0));
0057 
0058   toR_.clear();
0059   toZ_.clear();
0060 
0061   if (iSeed_ == 8 || iSeed_ == 9) {
0062     if (layer_ == 3) {
0063       rzmeanInv_[0] = 1.0 / settings_.rmean(2 - 1);
0064       rzmeanInv_[1] = 1.0 / settings_.rmean(3 - 1);
0065       rzmeanInv_[2] = 1.0 / settings_.rmean(4 - 1);
0066 
0067       rproj_[0] = settings_.rmean(0);
0068       rproj_[1] = settings_.rmean(4);
0069       rproj_[2] = settings_.rmean(5);
0070       lproj_[0] = 1;
0071       lproj_[1] = 5;
0072       lproj_[2] = 6;
0073 
0074       dproj_[0] = 1;
0075       dproj_[1] = 2;
0076       dproj_[2] = 0;
0077       toZ_.push_back(settings_.zmean(0));
0078       toZ_.push_back(settings_.zmean(1));
0079     }
0080     if (layer_ == 5) {
0081       rzmeanInv_[0] = 1.0 / settings_.rmean(4 - 1);
0082       rzmeanInv_[1] = 1.0 / settings_.rmean(5 - 1);
0083       rzmeanInv_[2] = 1.0 / settings_.rmean(6 - 1);
0084 
0085       rproj_[0] = settings_.rmean(0);
0086       rproj_[1] = settings_.rmean(1);
0087       rproj_[2] = settings_.rmean(2);
0088       lproj_[0] = 1;
0089       lproj_[1] = 2;
0090       lproj_[2] = 3;
0091 
0092       dproj_[0] = 0;
0093       dproj_[1] = 0;
0094       dproj_[2] = 0;
0095     }
0096     for (unsigned int i = 0; i < N_LAYER - 3; ++i)
0097       toR_.push_back(rproj_[i]);
0098   }
0099 
0100   if (iSeed_ == 10 || iSeed_ == 11) {
0101     if (layer_ == 2) {
0102       rzmeanInv_[0] = 1.0 / settings_.rmean(2 - 1);
0103       rzmeanInv_[1] = 1.0 / settings_.rmean(3 - 1);
0104       rzmeanInv_[2] = 1.0 / settings_.zmean(1 - 1);
0105 
0106       rproj_[0] = settings_.rmean(0);
0107       lproj_[0] = 1;
0108       lproj_[1] = -1;
0109       lproj_[2] = -1;
0110 
0111       zproj_[0] = settings_.zmean(1);
0112       zproj_[1] = settings_.zmean(2);
0113       zproj_[2] = settings_.zmean(3);
0114       dproj_[0] = 2;
0115       dproj_[1] = 3;
0116       dproj_[2] = 4;
0117     }
0118     if (disk_ == 1) {
0119       rzmeanInv_[0] = 1.0 / settings_.rmean(2 - 1);
0120       rzmeanInv_[1] = 1.0 / settings_.zmean(1 - 1);
0121       rzmeanInv_[2] = 1.0 / settings_.zmean(2 - 1);
0122 
0123       rproj_[0] = settings_.rmean(0);
0124       lproj_[0] = 1;
0125       lproj_[1] = -1;
0126       lproj_[2] = -1;
0127 
0128       zproj_[0] = settings_.zmean(2);
0129       zproj_[1] = settings_.zmean(3);
0130       zproj_[2] = settings_.zmean(4);
0131       dproj_[0] = 3;
0132       dproj_[1] = 4;
0133       dproj_[2] = 5;
0134     }
0135     toR_.push_back(settings_.rmean(0));
0136     for (unsigned int i = 0; i < N_DISK - 2; ++i)
0137       toZ_.push_back(zproj_[i]);
0138   }
0139 }
0140 
0141 void TrackletCalculatorDisplaced::addOutputProjection(TrackletProjectionsMemory*& outputProj, MemoryBase* memory) {
0142   outputProj = dynamic_cast<TrackletProjectionsMemory*>(memory);
0143   assert(outputProj != nullptr);
0144 }
0145 
0146 void TrackletCalculatorDisplaced::addOutput(MemoryBase* memory, string output) {
0147   if (settings_.writetrace()) {
0148     edm::LogVerbatim("Tracklet") << "In " << name_ << " adding output to " << memory->getName() << " to output "
0149                                  << output;
0150   }
0151 
0152   if (output == "trackpar") {
0153     auto* tmp = dynamic_cast<TrackletParametersMemory*>(memory);
0154     assert(tmp != nullptr);
0155     trackletpars_ = tmp;
0156     return;
0157   }
0158 
0159   if (output.substr(0, 7) == "projout") {
0160     //output is on the form 'projoutL2PHIC' or 'projoutD3PHIB'
0161     auto* tmp = dynamic_cast<TrackletProjectionsMemory*>(memory);
0162     assert(tmp != nullptr);
0163 
0164     unsigned int layerdisk = output[8] - '1';   //layer or disk counting from 0
0165     unsigned int phiregion = output[12] - 'A';  //phiregion counting from 0
0166 
0167     if (output[7] == 'L') {
0168       assert(layerdisk < N_LAYER);
0169       assert(phiregion < trackletprojlayers_[layerdisk].size());
0170       //check that phiregion not already initialized
0171       assert(trackletprojlayers_[layerdisk][phiregion] == nullptr);
0172       trackletprojlayers_[layerdisk][phiregion] = tmp;
0173       return;
0174     }
0175 
0176     if (output[7] == 'D') {
0177       assert(layerdisk < N_DISK);
0178       assert(phiregion < trackletprojdisks_[layerdisk].size());
0179       //check that phiregion not already initialized
0180       assert(trackletprojdisks_[layerdisk][phiregion] == nullptr);
0181       trackletprojdisks_[layerdisk][phiregion] = tmp;
0182       return;
0183     }
0184   }
0185 
0186   throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " Could not find output : " << output;
0187 }
0188 
0189 void TrackletCalculatorDisplaced::addInput(MemoryBase* memory, string input) {
0190   if (settings_.writetrace()) {
0191     edm::LogVerbatim("Tracklet") << "In " << name_ << " adding input from " << memory->getName() << " to input "
0192                                  << input;
0193   }
0194 
0195   if (input == "thirdallstubin") {
0196     auto* tmp = dynamic_cast<AllStubsMemory*>(memory);
0197     assert(tmp != nullptr);
0198     innerallstubs_.push_back(tmp);
0199     return;
0200   }
0201   if (input == "firstallstubin") {
0202     auto* tmp = dynamic_cast<AllStubsMemory*>(memory);
0203     assert(tmp != nullptr);
0204     middleallstubs_.push_back(tmp);
0205     return;
0206   }
0207   if (input == "secondallstubin") {
0208     auto* tmp = dynamic_cast<AllStubsMemory*>(memory);
0209     assert(tmp != nullptr);
0210     outerallstubs_.push_back(tmp);
0211     return;
0212   }
0213   if (input.find("stubtriplet") == 0) {
0214     auto* tmp = dynamic_cast<StubTripletsMemory*>(memory);
0215     assert(tmp != nullptr);
0216     stubtriplets_.push_back(tmp);
0217     return;
0218   }
0219   throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " Could not find input : " << input;
0220 }
0221 
0222 void TrackletCalculatorDisplaced::execute(unsigned int iSector, double phimin, double phimax) {
0223   unsigned int countall = 0;
0224   unsigned int countsel = 0;
0225 
0226   phimin_ = phimin;
0227   phimax_ = phimax;
0228   iSector_ = iSector;
0229 
0230   for (auto& stubtriplet : stubtriplets_) {
0231     if (trackletpars_->nTracklets() >= settings_.ntrackletmax()) {
0232       edm::LogVerbatim("Tracklet") << "Will break on too many tracklets in " << getName();
0233       break;
0234     }
0235     for (unsigned int i = 0; i < stubtriplet->nStubTriplets(); i++) {
0236       countall++;
0237 
0238       const Stub* innerFPGAStub = stubtriplet->getFPGAStub1(i);
0239       const L1TStub* innerStub = innerFPGAStub->l1tstub();
0240 
0241       const Stub* middleFPGAStub = stubtriplet->getFPGAStub2(i);
0242       const L1TStub* middleStub = middleFPGAStub->l1tstub();
0243 
0244       const Stub* outerFPGAStub = stubtriplet->getFPGAStub3(i);
0245       const L1TStub* outerStub = outerFPGAStub->l1tstub();
0246 
0247       if (settings_.debugTracklet())
0248         edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced execute " << getName() << "[" << iSector_ << "]";
0249 
0250       if (innerFPGAStub->layerdisk() < N_LAYER && middleFPGAStub->layerdisk() < N_LAYER &&
0251           outerFPGAStub->layerdisk() < N_LAYER) {
0252         //barrel+barrel seeding
0253         bool accept = LLLSeeding(innerFPGAStub, innerStub, middleFPGAStub, middleStub, outerFPGAStub, outerStub);
0254         if (accept)
0255           countsel++;
0256       } else if (innerFPGAStub->layerdisk() >= N_LAYER && middleFPGAStub->layerdisk() >= N_LAYER &&
0257                  outerFPGAStub->layerdisk() >= N_LAYER) {
0258         throw cms::Exception("LogicError") << __FILE__ << " " << __LINE__ << " Invalid seeding!";
0259       } else {
0260         //layer+disk seeding
0261         if (innerFPGAStub->layerdisk() < N_LAYER && middleFPGAStub->layerdisk() >= N_LAYER &&
0262             outerFPGAStub->layerdisk() >= N_LAYER) {  //D1D2L2
0263           bool accept = DDLSeeding(innerFPGAStub, innerStub, middleFPGAStub, middleStub, outerFPGAStub, outerStub);
0264           if (accept)
0265             countsel++;
0266         } else if (innerFPGAStub->layerdisk() >= N_LAYER && middleFPGAStub->layerdisk() < N_LAYER &&
0267                    outerFPGAStub->layerdisk() < N_LAYER) {  //L2L3D1
0268           bool accept = LLDSeeding(innerFPGAStub, innerStub, middleFPGAStub, middleStub, outerFPGAStub, outerStub);
0269           if (accept)
0270             countsel++;
0271         } else {
0272           throw cms::Exception("LogicError") << __FILE__ << " " << __LINE__ << " Invalid seeding!";
0273         }
0274       }
0275 
0276       if (trackletpars_->nTracklets() >= settings_.ntrackletmax()) {
0277         edm::LogVerbatim("Tracklet") << "Will break on number of tracklets in " << getName();
0278         break;
0279       }
0280 
0281       if (countall >= settings_.maxStep("TC")) {
0282         if (settings_.debugTracklet())
0283           edm::LogVerbatim("Tracklet") << "Will break on MAXTC 1";
0284         break;
0285       }
0286       if (settings_.debugTracklet())
0287         edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced execute done";
0288     }
0289     if (countall >= settings_.maxStep("TC")) {
0290       if (settings_.debugTracklet())
0291         edm::LogVerbatim("Tracklet") << "Will break on MAXTC 2";
0292       break;
0293     }
0294   }
0295 
0296   if (settings_.writeMonitorData("TPD")) {
0297     globals_->ofstream("trackletcalculatordisplaced.txt") << getName() << " " << countall << " " << countsel << endl;
0298   }
0299 }
0300 
0301 void TrackletCalculatorDisplaced::addDiskProj(Tracklet* tracklet, int disk) {
0302   disk = std::abs(disk);
0303   FPGAWord fpgar = tracklet->proj(N_LAYER + disk - 1).fpgarzproj();
0304 
0305   if (fpgar.value() * settings_.krprojshiftdisk() < settings_.rmindiskvm())
0306     return;
0307   if (fpgar.value() * settings_.krprojshiftdisk() > settings_.rmaxdisk())
0308     return;
0309 
0310   FPGAWord fpgaphi = tracklet->proj(N_LAYER + disk - 1).fpgaphiproj();
0311 
0312   int iphivmRaw = fpgaphi.value() >> (fpgaphi.nbits() - 5);
0313   int iphi = iphivmRaw / (32 / settings_.nallstubs(disk + N_LAYER - 1));
0314 
0315   addProjectionDisk(disk, iphi, trackletprojdisks_[disk - 1][iphi], tracklet);
0316 }
0317 
0318 bool TrackletCalculatorDisplaced::addLayerProj(Tracklet* tracklet, int layer) {
0319   assert(layer > 0);
0320 
0321   FPGAWord fpgaz = tracklet->proj(layer - 1).fpgarzproj();
0322   FPGAWord fpgaphi = tracklet->proj(layer - 1).fpgaphiproj();
0323 
0324   if (fpgaz.atExtreme())
0325     return false;
0326 
0327   if (std::abs(fpgaz.value() * settings_.kz()) > settings_.zlength())
0328     return false;
0329 
0330   int iphivmRaw = fpgaphi.value() >> (fpgaphi.nbits() - 5);
0331   int iphi = iphivmRaw / (32 / settings_.nallstubs(layer - 1));
0332 
0333   addProjection(layer, iphi, trackletprojlayers_[layer - 1][iphi], tracklet);
0334 
0335   return true;
0336 }
0337 
0338 void TrackletCalculatorDisplaced::addProjection(int layer,
0339                                                 int iphi,
0340                                                 TrackletProjectionsMemory* trackletprojs,
0341                                                 Tracklet* tracklet) {
0342   if (trackletprojs == nullptr) {
0343     if (settings_.warnNoMem()) {
0344       edm::LogVerbatim("Tracklet") << "No projection memory exists in " << getName() << " for layer = " << layer
0345                                    << " iphi = " << iphi + 1;
0346     }
0347     return;
0348   }
0349   assert(trackletprojs != nullptr);
0350   trackletprojs->addProj(tracklet);
0351 }
0352 
0353 void TrackletCalculatorDisplaced::addProjectionDisk(int disk,
0354                                                     int iphi,
0355                                                     TrackletProjectionsMemory* trackletprojs,
0356                                                     Tracklet* tracklet) {
0357   if (trackletprojs == nullptr) {
0358     if (layer_ == 3 && abs(disk) == 3)
0359       return;  //L3L4 projections to D3 are not used.
0360     if (settings_.warnNoMem()) {
0361       edm::LogVerbatim("Tracklet") << "No projection memory exists in " << getName() << " for disk = " << abs(disk)
0362                                    << " iphi = " << iphi + 1;
0363     }
0364     return;
0365   }
0366   assert(trackletprojs != nullptr);
0367   if (settings_.debugTracklet())
0368     edm::LogVerbatim("Tracklet") << getName() << " adding projection to " << trackletprojs->getName();
0369   trackletprojs->addProj(tracklet);
0370 }
0371 
0372 bool TrackletCalculatorDisplaced::LLLSeeding(const Stub* innerFPGAStub,
0373                                              const L1TStub* innerStub,
0374                                              const Stub* middleFPGAStub,
0375                                              const L1TStub* middleStub,
0376                                              const Stub* outerFPGAStub,
0377                                              const L1TStub* outerStub) {
0378   if (settings_.debugTracklet())
0379     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName() << " " << layer_
0380                                  << " trying stub triplet in layer (L L L): " << innerFPGAStub->layer().value() << " "
0381                                  << middleFPGAStub->layer().value() << " " << outerFPGAStub->layer().value();
0382 
0383   assert(outerFPGAStub->layerdisk() < N_LAYER);
0384 
0385   double r1 = innerStub->r();
0386   double z1 = innerStub->z();
0387   double phi1 = innerStub->phi();
0388 
0389   double r2 = middleStub->r();
0390   double z2 = middleStub->z();
0391   double phi2 = middleStub->phi();
0392 
0393   double r3 = outerStub->r();
0394   double z3 = outerStub->z();
0395   double phi3 = outerStub->phi();
0396 
0397   int take3 = 0;
0398   if (layer_ == 5)
0399     take3 = 1;
0400   unsigned ndisks = 0;
0401 
0402   double rinv, phi0, d0, t, z0;
0403 
0404   Projection projs[N_LAYER + N_DISK];
0405 
0406   double phiproj[N_LAYER - 2], zproj[N_LAYER - 2], phider[N_LAYER - 2], zder[N_LAYER - 2];
0407   double phiprojdisk[N_DISK], rprojdisk[N_DISK], phiderdisk[N_DISK], rderdisk[N_DISK];
0408 
0409   exacttracklet(r1,
0410                 z1,
0411                 phi1,
0412                 r2,
0413                 z2,
0414                 phi2,
0415                 r3,
0416                 z3,
0417                 phi3,
0418                 take3,
0419                 rinv,
0420                 phi0,
0421                 d0,
0422                 t,
0423                 z0,
0424                 phiproj,
0425                 zproj,
0426                 phiprojdisk,
0427                 rprojdisk,
0428                 phider,
0429                 zder,
0430                 phiderdisk,
0431                 rderdisk);
0432   if (settings_.debugTracklet())
0433     edm::LogVerbatim("Tracklet") << __LINE__ << ":" << __FILE__ << "LLL Exact values " << innerFPGAStub->isBarrel()
0434                                  << middleFPGAStub->isBarrel() << outerFPGAStub->isBarrel() << " " << phi1 << ", " << z1
0435                                  << ", " << r1 << ", " << phi2 << ", " << z2 << ", " << r2 << ", " << phi3 << ", " << z3
0436                                  << ", " << r3 << endl;
0437 
0438   if (settings_.useapprox()) {
0439     phi1 = innerFPGAStub->phiapprox(phimin_, phimax_);
0440     z1 = innerFPGAStub->zapprox();
0441     r1 = innerFPGAStub->rapprox();
0442 
0443     phi2 = middleFPGAStub->phiapprox(phimin_, phimax_);
0444     z2 = middleFPGAStub->zapprox();
0445     r2 = middleFPGAStub->rapprox();
0446 
0447     phi3 = outerFPGAStub->phiapprox(phimin_, phimax_);
0448     z3 = outerFPGAStub->zapprox();
0449     r3 = outerFPGAStub->rapprox();
0450   }
0451 
0452   if (settings_.debugTracklet())
0453     edm::LogVerbatim("Tracklet") << __LINE__ << ":" << __FILE__ << "LLL Approx values " << innerFPGAStub->isBarrel()
0454                                  << middleFPGAStub->isBarrel() << outerFPGAStub->isBarrel() << " " << phi1 << ", " << z1
0455                                  << ", " << r1 << ", " << phi2 << ", " << z2 << ", " << r2 << ", " << phi3 << ", " << z3
0456                                  << ", " << r3 << endl;
0457 
0458   double rinvapprox, phi0approx, d0approx, tapprox, z0approx;
0459   double phiprojapprox[N_LAYER - 2], zprojapprox[N_LAYER - 2], phiderapprox[N_LAYER - 2], zderapprox[N_LAYER - 2];
0460   double phiprojdiskapprox[N_DISK], rprojdiskapprox[N_DISK];
0461   double phiderdiskapprox[N_DISK], rderdiskapprox[N_DISK];
0462 
0463   //TODO: implement the actual integer calculation
0464   if (settings_.useapprox()) {
0465     approxtracklet(r1,
0466                    z1,
0467                    phi1,
0468                    r2,
0469                    z2,
0470                    phi2,
0471                    r3,
0472                    z3,
0473                    phi3,
0474                    take3,
0475                    ndisks,
0476                    rinvapprox,
0477                    phi0approx,
0478                    d0approx,
0479                    tapprox,
0480                    z0approx,
0481                    phiprojapprox,
0482                    zprojapprox,
0483                    phiderapprox,
0484                    zderapprox,
0485                    phiprojdiskapprox,
0486                    rprojdiskapprox,
0487                    phiderdiskapprox,
0488                    rderdiskapprox);
0489   } else {
0490     rinvapprox = rinv;
0491     phi0approx = phi0;
0492     d0approx = d0;
0493     tapprox = t;
0494     z0approx = z0;
0495 
0496     for (unsigned int i = 0; i < toR_.size(); ++i) {
0497       phiprojapprox[i] = phiproj[i];
0498       zprojapprox[i] = zproj[i];
0499       phiderapprox[i] = phider[i];
0500       zderapprox[i] = zder[i];
0501     }
0502 
0503     for (unsigned int i = 0; i < toZ_.size(); ++i) {
0504       phiprojdiskapprox[i] = phiprojdisk[i];
0505       rprojdiskapprox[i] = rprojdisk[i];
0506       phiderdiskapprox[i] = phiderdisk[i];
0507       rderdiskapprox[i] = rderdisk[i];
0508     }
0509   }
0510 
0511   //store the approcximate results
0512 
0513   if (settings_.debugTracklet()) {
0514     edm::LogVerbatim("Tracklet") << "rinvapprox: " << rinvapprox << " rinv: " << rinv << endl;
0515     edm::LogVerbatim("Tracklet") << "phi0approx: " << phi0approx << " phi0: " << phi0 << endl;
0516     edm::LogVerbatim("Tracklet") << "d0approx: " << d0approx << " d0: " << d0 << endl;
0517     edm::LogVerbatim("Tracklet") << "tapprox: " << tapprox << " t: " << t << endl;
0518     edm::LogVerbatim("Tracklet") << "z0approx: " << z0approx << " z0: " << z0 << endl;
0519   }
0520 
0521   for (unsigned int i = 0; i < toR_.size(); ++i) {
0522     if (settings_.debugTracklet()) {
0523       edm::LogVerbatim("Tracklet") << "phiprojapprox[" << i << "]: " << phiprojapprox[i] << " phiproj[" << i
0524                                    << "]: " << phiproj[i] << endl;
0525       edm::LogVerbatim("Tracklet") << "zprojapprox[" << i << "]: " << zprojapprox[i] << " zproj[" << i
0526                                    << "]: " << zproj[i] << endl;
0527       edm::LogVerbatim("Tracklet") << "phiderapprox[" << i << "]: " << phiderapprox[i] << " phider[" << i
0528                                    << "]: " << phider[i] << endl;
0529       edm::LogVerbatim("Tracklet") << "zderapprox[" << i << "]: " << zderapprox[i] << " zder[" << i << "]: " << zder[i]
0530                                    << endl;
0531     }
0532   }
0533 
0534   for (unsigned int i = 0; i < toZ_.size(); ++i) {
0535     if (settings_.debugTracklet()) {
0536       edm::LogVerbatim("Tracklet") << "phiprojdiskapprox[" << i << "]: " << phiprojdiskapprox[i] << " phiprojdisk[" << i
0537                                    << "]: " << phiprojdisk[i] << endl;
0538       edm::LogVerbatim("Tracklet") << "rprojdiskapprox[" << i << "]: " << rprojdiskapprox[i] << " rprojdisk[" << i
0539                                    << "]: " << rprojdisk[i] << endl;
0540       edm::LogVerbatim("Tracklet") << "phiderdiskapprox[" << i << "]: " << phiderdiskapprox[i] << " phiderdisk[" << i
0541                                    << "]: " << phiderdisk[i] << endl;
0542       edm::LogVerbatim("Tracklet") << "rderdiskapprox[" << i << "]: " << rderdiskapprox[i] << " rderdisk[" << i
0543                                    << "]: " << rderdisk[i] << endl;
0544     }
0545   }
0546 
0547   //now binary
0548   double krinv = settings_.kphi1() / settings_.kr() * pow(2, settings_.rinv_shift()),
0549          kphi0 = settings_.kphi1() * pow(2, settings_.phi0_shift()),
0550          kt = settings_.kz() / settings_.kr() * pow(2, settings_.t_shift()),
0551          kz0 = settings_.kz() * pow(2, settings_.z0_shift()),
0552          kphiproj = settings_.kphi1() * pow(2, settings_.SS_phiL_shift()),
0553          kphider = settings_.kphi1() / settings_.kr() * pow(2, settings_.SS_phiderL_shift()),
0554          kzproj = settings_.kz() * pow(2, settings_.PS_zL_shift()),
0555          kzder = settings_.kz() / settings_.kr() * pow(2, settings_.PS_zderL_shift()),
0556          kphiprojdisk = settings_.kphi1() * pow(2, settings_.SS_phiD_shift()),
0557          kphiderdisk = settings_.kphi1() / settings_.kr() * pow(2, settings_.SS_phiderD_shift()),
0558          krprojdisk = settings_.kr() * pow(2, settings_.PS_rD_shift()),
0559          krderdisk = settings_.kr() / settings_.kz() * pow(2, settings_.PS_rderD_shift());
0560 
0561   int irinv, iphi0, id0, it, iz0;
0562   int iphiproj[N_LAYER - 2], izproj[N_LAYER - 2], iphider[N_LAYER - 2], izder[N_LAYER - 2];
0563   int iphiprojdisk[N_DISK], irprojdisk[N_DISK], iphiderdisk[N_DISK], irderdisk[N_DISK];
0564 
0565   //store the binary results
0566   irinv = rinvapprox / krinv;
0567   iphi0 = phi0approx / kphi0;
0568   id0 = d0approx / settings_.kd0();
0569   it = tapprox / kt;
0570   iz0 = z0approx / kz0;
0571 
0572   bool success = true;
0573   if (std::abs(rinvapprox) > settings_.rinvcut()) {
0574     if (settings_.debugTracklet())
0575       edm::LogVerbatim("Tracklet") << "TrackletCalculator::LLL Seeding irinv too large: " << rinvapprox << "(" << irinv
0576                                    << ")";
0577     success = false;
0578   }
0579   if (std::abs(z0approx) > settings_.disp_z0cut()) {
0580     if (settings_.debugTracklet())
0581       edm::LogVerbatim("Tracklet") << "Failed tracklet z0 cut " << z0approx << " in layer " << layer_;
0582     success = false;
0583   }
0584   if (std::abs(d0approx) > settings_.maxd0()) {
0585     if (settings_.debugTracklet())
0586       edm::LogVerbatim("Tracklet") << "Failed tracklet approx d0 cut " << d0approx;
0587     success = false;
0588   }
0589   if (std::abs(d0) > settings_.maxd0()) {
0590     if (settings_.debugTracklet())
0591       edm::LogVerbatim("Tracklet") << "Failed tracklet exact d0 cut " << d0;
0592     success = false;
0593   }
0594 
0595   if (!success) {
0596     return false;
0597   }
0598 
0599   double phicritapprox = phi0approx - asin((0.5 * settings_.rcrit() * rinvapprox) + (d0approx / settings_.rcrit()));
0600   int phicrit = iphi0 - 2 * irinv - 2 * id0;
0601 
0602   int iphicritmincut = settings_.phicritminmc() / globals_->ITC_L1L2()->phi0_final.K();
0603   int iphicritmaxcut = settings_.phicritmaxmc() / globals_->ITC_L1L2()->phi0_final.K();
0604 
0605   bool keepapprox = (phicritapprox > settings_.phicritminmc()) && (phicritapprox < settings_.phicritmaxmc()),
0606        keep = (phicrit > iphicritmincut) && (phicrit < iphicritmaxcut);
0607 
0608   if (settings_.debugTracklet())
0609     if (keep && !keepapprox)
0610       edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced::LLLSeeding tracklet kept with exact phicrit cut "
0611                                       "but not approximate, phicritapprox: "
0612                                    << phicritapprox;
0613   if (settings_.usephicritapprox()) {
0614     if (!keepapprox) {
0615       return false;
0616     }
0617   } else {
0618     if (!keep) {
0619       return false;
0620     }
0621   }
0622 
0623   for (unsigned int i = 0; i < toR_.size(); ++i) {
0624     iphiproj[i] = phiprojapprox[i] / kphiproj;
0625     izproj[i] = zprojapprox[i] / kzproj;
0626 
0627     iphider[i] = phiderapprox[i] / kphider;
0628     izder[i] = zderapprox[i] / kzder;
0629 
0630     //check that z projection is in range
0631     if (izproj[i] < -(1 << (settings_.nzbitsstub(0) - 1)))
0632       continue;
0633     if (izproj[i] >= (1 << (settings_.nzbitsstub(0) - 1)))
0634       continue;
0635 
0636     //check that phi projection is in range
0637     if (iphiproj[i] >= (1 << settings_.nphibitsstub(N_LAYER - 1)) - 1)
0638       continue;
0639     if (iphiproj[i] <= 0)
0640       continue;
0641 
0642     //adjust number of bits for phi and z projection
0643     if (rproj_[i] < settings_.rPS2S()) {
0644       iphiproj[i] >>= (settings_.nphibitsstub(N_LAYER - 1) - settings_.nphibitsstub(0));
0645       if (iphiproj[i] >= (1 << settings_.nphibitsstub(0)) - 1)
0646         iphiproj[i] = (1 << settings_.nphibitsstub(0)) - 2;  //-2 not to hit atExtreme
0647     } else {
0648       izproj[i] >>= (settings_.nzbitsstub(0) - settings_.nzbitsstub(N_LAYER - 1));
0649     }
0650 
0651     if (rproj_[i] < settings_.rPS2S()) {
0652       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL123() - 1))) {
0653         iphider[i] = -(1 << (settings_.nbitsphiprojderL123() - 1));
0654       }
0655       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL123() - 1))) {
0656         iphider[i] = (1 << (settings_.nbitsphiprojderL123() - 1)) - 1;
0657       }
0658     } else {
0659       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL456() - 1))) {
0660         iphider[i] = -(1 << (settings_.nbitsphiprojderL456() - 1));
0661       }
0662       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL456() - 1))) {
0663         iphider[i] = (1 << (settings_.nbitsphiprojderL456() - 1)) - 1;
0664       }
0665     }
0666 
0667     projs[lproj_[i] - 1].init(settings_,
0668                               lproj_[i] - 1,
0669                               iphiproj[i],
0670                               izproj[i],
0671                               iphider[i],
0672                               izder[i],
0673                               phiproj[i],
0674                               zproj[i],
0675                               phider[i],
0676                               zder[i],
0677                               phiprojapprox[i],
0678                               zprojapprox[i],
0679                               phiderapprox[i],
0680                               zderapprox[i],
0681                               false);
0682   }
0683 
0684   if (std::abs(it * kt) > 1.0) {
0685     for (unsigned int i = 0; i < toZ_.size(); ++i) {
0686       iphiprojdisk[i] = phiprojdiskapprox[i] / kphiprojdisk;
0687       irprojdisk[i] = rprojdiskapprox[i] / krprojdisk;
0688 
0689       iphiderdisk[i] = phiderdiskapprox[i] / kphiderdisk;
0690       irderdisk[i] = rderdiskapprox[i] / krderdisk;
0691 
0692       //check phi projection in range
0693       if (iphiprojdisk[i] <= 0)
0694         continue;
0695       if (iphiprojdisk[i] >= (1 << settings_.nphibitsstub(0)) - 1)
0696         continue;
0697 
0698       //check r projection in range
0699       if (rprojdiskapprox[i] < settings_.rmindisk() || rprojdiskapprox[i] >= settings_.rmaxdisk())
0700         continue;
0701 
0702       projs[N_LAYER + i].init(settings_,
0703                               N_LAYER + i,
0704                               iphiprojdisk[i],
0705                               irprojdisk[i],
0706                               iphiderdisk[i],
0707                               irderdisk[i],
0708                               phiprojdisk[i],
0709                               rprojdisk[i],
0710                               phiderdisk[i],
0711                               rderdisk[i],
0712                               phiprojdiskapprox[i],
0713                               rprojdiskapprox[i],
0714                               phiderdisk[i],
0715                               rderdisk[i],
0716                               false);
0717     }
0718   }
0719 
0720   if (settings_.writeMonitorData("TrackletPars")) {
0721     globals_->ofstream("trackletpars.txt")
0722         << layer_ << " , " << rinv << " , " << rinvapprox << " , " << phi0 << " , " << phi0approx << " , " << t << " , "
0723         << tapprox << " , " << z0 << " , " << z0approx << " , " << d0 << " , " << d0approx << endl;
0724   }
0725 
0726   Tracklet* tracklet = new Tracklet(settings_,
0727                                     iSeed_,
0728                                     innerFPGAStub,
0729                                     middleFPGAStub,
0730                                     outerFPGAStub,
0731                                     rinv,
0732                                     phi0,
0733                                     d0,
0734                                     z0,
0735                                     t,
0736                                     rinvapprox,
0737                                     phi0approx,
0738                                     d0approx,
0739                                     z0approx,
0740                                     tapprox,
0741                                     irinv,
0742                                     iphi0,
0743                                     id0,
0744                                     iz0,
0745                                     it,
0746                                     projs,
0747                                     false);
0748 
0749   if (settings_.debugTracklet())
0750     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName()
0751                                  << " Found LLL tracklet in sector = " << iSector_ << " phi0 = " << phi0;
0752 
0753   tracklet->setTrackletIndex(trackletpars_->nTracklets());
0754   tracklet->setTCIndex(TCIndex_);
0755 
0756   if (settings_.writeMonitorData("Seeds")) {
0757     ofstream fout("seeds.txt", ofstream::app);
0758     fout << __FILE__ << ":" << __LINE__ << " " << name_ << "_" << iSector_ << " " << tracklet->getISeed() << endl;
0759     fout.close();
0760   }
0761   trackletpars_->addTracklet(tracklet);
0762 
0763   bool addL5 = false;
0764   bool addL6 = false;
0765   for (unsigned int j = 0; j < toR_.size(); j++) {
0766     bool added = false;
0767 
0768     if (settings_.debugTracklet())
0769       edm::LogVerbatim("Tracklet") << "adding layer projection " << j << "/" << toR_.size() << " " << lproj_[j];
0770     if (tracklet->validProj(lproj_[j] - 1)) {
0771       added = addLayerProj(tracklet, lproj_[j]);
0772       if (added && lproj_[j] == 5)
0773         addL5 = true;
0774       if (added && lproj_[j] == 6)
0775         addL6 = true;
0776     }
0777   }
0778 
0779   for (unsigned int j = 0; j < toZ_.size(); j++) {
0780     int disk = dproj_[j];
0781     if (disk == 0)
0782       continue;
0783     if (disk == 2 && addL5)
0784       continue;
0785     if (disk == 1 && addL6)
0786       continue;
0787     if (it < 0)
0788       disk = -disk;
0789     if (settings_.debugTracklet())
0790       edm::LogVerbatim("Tracklet") << "adding disk projection " << j << "/" << toZ_.size() << " " << disk;
0791     if (tracklet->validProj(N_LAYER + abs(disk) - 1)) {
0792       addDiskProj(tracklet, disk);
0793     }
0794   }
0795 
0796   return true;
0797 }
0798 
0799 bool TrackletCalculatorDisplaced::DDLSeeding(const Stub* innerFPGAStub,
0800                                              const L1TStub* innerStub,
0801                                              const Stub* middleFPGAStub,
0802                                              const L1TStub* middleStub,
0803                                              const Stub* outerFPGAStub,
0804                                              const L1TStub* outerStub) {
0805   if (settings_.debugTracklet())
0806     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName() << " " << layer_
0807                                  << " trying stub triplet in  (L2 D1 D2): " << innerFPGAStub->layer().value() << " "
0808                                  << middleFPGAStub->disk().value() << " " << outerFPGAStub->disk().value();
0809 
0810   int take3 = 1;  //D1D2L2
0811   unsigned ndisks = 2;
0812 
0813   double r1 = innerStub->r();
0814   double z1 = innerStub->z();
0815   double phi1 = innerStub->phi();
0816 
0817   double r2 = middleStub->r();
0818   double z2 = middleStub->z();
0819   double phi2 = middleStub->phi();
0820 
0821   double r3 = outerStub->r();
0822   double z3 = outerStub->z();
0823   double phi3 = outerStub->phi();
0824 
0825   double rinv, phi0, d0, t, z0;
0826 
0827   double phiproj[N_LAYER - 2], zproj[N_LAYER - 2], phider[N_LAYER - 2], zder[N_LAYER - 2];
0828   double phiprojdisk[N_DISK], rprojdisk[N_DISK], phiderdisk[N_DISK], rderdisk[N_DISK];
0829 
0830   exacttracklet(r1,
0831                 z1,
0832                 phi1,
0833                 r2,
0834                 z2,
0835                 phi2,
0836                 r3,
0837                 z3,
0838                 phi3,
0839                 take3,
0840                 rinv,
0841                 phi0,
0842                 d0,
0843                 t,
0844                 z0,
0845                 phiproj,
0846                 zproj,
0847                 phiprojdisk,
0848                 rprojdisk,
0849                 phider,
0850                 zder,
0851                 phiderdisk,
0852                 rderdisk);
0853 
0854   if (settings_.debugTracklet())
0855     edm::LogVerbatim("Tracklet") << __LINE__ << ":" << __FILE__ << " DLL Exact values " << innerFPGAStub->isBarrel()
0856                                  << middleFPGAStub->isBarrel() << outerFPGAStub->isBarrel() << " " << phi1 << ", " << z1
0857                                  << ", " << r1 << ", " << phi2 << ", " << z2 << ", " << r2 << ", " << phi3 << ", " << z3
0858                                  << ", " << r3 << endl;
0859 
0860   if (settings_.useapprox()) {
0861     phi1 = innerFPGAStub->phiapprox(phimin_, phimax_);
0862     z1 = innerFPGAStub->zapprox();
0863     r1 = innerFPGAStub->rapprox();
0864 
0865     phi2 = middleFPGAStub->phiapprox(phimin_, phimax_);
0866     z2 = middleFPGAStub->zapprox();
0867     r2 = middleFPGAStub->rapprox();
0868 
0869     phi3 = outerFPGAStub->phiapprox(phimin_, phimax_);
0870     z3 = outerFPGAStub->zapprox();
0871     r3 = outerFPGAStub->rapprox();
0872   }
0873 
0874   if (settings_.debugTracklet())
0875     edm::LogVerbatim("Tracklet") << __LINE__ << ":" << __FILE__ << "DLL Approx values " << innerFPGAStub->isBarrel()
0876                                  << middleFPGAStub->isBarrel() << outerFPGAStub->isBarrel() << " " << phi1 << ", " << z1
0877                                  << ", " << r1 << ", " << phi2 << ", " << z2 << ", " << r2 << ", " << phi3 << ", " << z3
0878                                  << ", " << r3 << endl;
0879 
0880   double rinvapprox, phi0approx, d0approx, tapprox, z0approx;
0881   double phiprojapprox[N_LAYER - 2], zprojapprox[N_LAYER - 2], phiderapprox[N_LAYER - 2], zderapprox[N_LAYER - 2];
0882   double phiprojdiskapprox[N_DISK], rprojdiskapprox[N_DISK];
0883   double phiderdiskapprox[N_DISK], rderdiskapprox[N_DISK];
0884 
0885   //TODO: implement the actual integer calculation
0886   if (settings_.useapprox()) {
0887     approxtracklet(r1,
0888                    z1,
0889                    phi1,
0890                    r2,
0891                    z2,
0892                    phi2,
0893                    r3,
0894                    z3,
0895                    phi3,
0896                    take3,
0897                    ndisks,
0898                    rinvapprox,
0899                    phi0approx,
0900                    d0approx,
0901                    tapprox,
0902                    z0approx,
0903                    phiprojapprox,
0904                    zprojapprox,
0905                    phiderapprox,
0906                    zderapprox,
0907                    phiprojdiskapprox,
0908                    rprojdiskapprox,
0909                    phiderdiskapprox,
0910                    rderdiskapprox);
0911   } else {
0912     rinvapprox = rinv;
0913     phi0approx = phi0;
0914     d0approx = d0;
0915     tapprox = t;
0916     z0approx = z0;
0917 
0918     for (unsigned int i = 0; i < toR_.size(); ++i) {
0919       phiprojapprox[i] = phiproj[i];
0920       zprojapprox[i] = zproj[i];
0921       phiderapprox[i] = phider[i];
0922       zderapprox[i] = zder[i];
0923     }
0924 
0925     for (unsigned int i = 0; i < toZ_.size(); ++i) {
0926       phiprojdiskapprox[i] = phiprojdisk[i];
0927       rprojdiskapprox[i] = rprojdisk[i];
0928       phiderdiskapprox[i] = phiderdisk[i];
0929       rderdiskapprox[i] = rderdisk[i];
0930     }
0931   }
0932 
0933   //store the approcximate results
0934   if (settings_.debugTracklet()) {
0935     edm::LogVerbatim("Tracklet") << "rinvapprox: " << rinvapprox << " rinv: " << rinv << endl;
0936     edm::LogVerbatim("Tracklet") << "phi0approx: " << phi0approx << " phi0: " << phi0 << endl;
0937     edm::LogVerbatim("Tracklet") << "d0approx: " << d0approx << " d0: " << d0 << endl;
0938     edm::LogVerbatim("Tracklet") << "tapprox: " << tapprox << " t: " << t << endl;
0939     edm::LogVerbatim("Tracklet") << "z0approx: " << z0approx << " z0: " << z0 << endl;
0940   }
0941 
0942   for (unsigned int i = 0; i < toR_.size(); ++i) {
0943     if (settings_.debugTracklet()) {
0944       edm::LogVerbatim("Tracklet") << "phiprojapprox[" << i << "]: " << phiprojapprox[i] << " phiproj[" << i
0945                                    << "]: " << phiproj[i] << endl;
0946       edm::LogVerbatim("Tracklet") << "zprojapprox[" << i << "]: " << zprojapprox[i] << " zproj[" << i
0947                                    << "]: " << zproj[i] << endl;
0948       edm::LogVerbatim("Tracklet") << "phiderapprox[" << i << "]: " << phiderapprox[i] << " phider[" << i
0949                                    << "]: " << phider[i] << endl;
0950       edm::LogVerbatim("Tracklet") << "zderapprox[" << i << "]: " << zderapprox[i] << " zder[" << i << "]: " << zder[i]
0951                                    << endl;
0952     }
0953   }
0954 
0955   for (unsigned int i = 0; i < toZ_.size(); ++i) {
0956     if (settings_.debugTracklet()) {
0957       edm::LogVerbatim("Tracklet") << "phiprojdiskapprox[" << i << "]: " << phiprojdiskapprox[i] << " phiprojdisk[" << i
0958                                    << "]: " << phiprojdisk[i] << endl;
0959       edm::LogVerbatim("Tracklet") << "rprojdiskapprox[" << i << "]: " << rprojdiskapprox[i] << " rprojdisk[" << i
0960                                    << "]: " << rprojdisk[i] << endl;
0961       edm::LogVerbatim("Tracklet") << "phiderdiskapprox[" << i << "]: " << phiderdiskapprox[i] << " phiderdisk[" << i
0962                                    << "]: " << phiderdisk[i] << endl;
0963       edm::LogVerbatim("Tracklet") << "rderdiskapprox[" << i << "]: " << rderdiskapprox[i] << " rderdisk[" << i
0964                                    << "]: " << rderdisk[i] << endl;
0965     }
0966   }
0967 
0968   //now binary
0969   double krinv = settings_.kphi1() / settings_.kr() * pow(2, settings_.rinv_shift()),
0970          kphi0 = settings_.kphi1() * pow(2, settings_.phi0_shift()),
0971          kt = settings_.kz() / settings_.kr() * pow(2, settings_.t_shift()),
0972          kz0 = settings_.kz() * pow(2, settings_.z0_shift()),
0973          kphiproj = settings_.kphi1() * pow(2, settings_.SS_phiL_shift()),
0974          kphider = settings_.kphi1() / settings_.kr() * pow(2, settings_.SS_phiderL_shift()),
0975          kzproj = settings_.kz() * pow(2, settings_.PS_zL_shift()),
0976          kzder = settings_.kz() / settings_.kr() * pow(2, settings_.PS_zderL_shift()),
0977          kphiprojdisk = settings_.kphi1() * pow(2, settings_.SS_phiD_shift()),
0978          kphiderdisk = settings_.kphi1() / settings_.kr() * pow(2, settings_.SS_phiderD_shift()),
0979          krprojdisk = settings_.kr() * pow(2, settings_.PS_rD_shift()),
0980          krderdisk = settings_.kr() / settings_.kz() * pow(2, settings_.PS_rderD_shift());
0981 
0982   int irinv, iphi0, id0, it, iz0;
0983   int iphiproj[N_LAYER - 2], izproj[N_LAYER - 2], iphider[N_LAYER - 2], izder[N_LAYER - 2];
0984   int iphiprojdisk[N_DISK], irprojdisk[N_DISK], iphiderdisk[N_DISK], irderdisk[N_DISK];
0985 
0986   //store the binary results
0987   irinv = rinvapprox / krinv;
0988   iphi0 = phi0approx / kphi0;
0989   id0 = d0approx / settings_.kd0();
0990   it = tapprox / kt;
0991   iz0 = z0approx / kz0;
0992 
0993   bool success = true;
0994   if (std::abs(rinvapprox) > settings_.rinvcut()) {
0995     if (settings_.debugTracklet())
0996       edm::LogVerbatim("Tracklet") << "TrackletCalculator::DDL Seeding irinv too large: " << rinvapprox << "(" << irinv
0997                                    << ")";
0998     success = false;
0999   }
1000   if (std::abs(z0approx) > settings_.disp_z0cut()) {
1001     if (settings_.debugTracklet())
1002       edm::LogVerbatim("Tracklet") << "Failed tracklet z0 cut " << z0approx;
1003     success = false;
1004   }
1005   if (std::abs(d0approx) > settings_.maxd0()) {
1006     if (settings_.debugTracklet())
1007       edm::LogVerbatim("Tracklet") << "Failed tracklet approx d0 cut " << d0approx;
1008     success = false;
1009   }
1010   if (std::abs(d0) > settings_.maxd0()) {
1011     if (settings_.debugTracklet())
1012       edm::LogVerbatim("Tracklet") << "Failed tracklet exact d0 cut " << d0;
1013     success = false;
1014   }
1015 
1016   if (!success)
1017     return false;
1018 
1019   double phicritapprox = phi0approx - asin((0.5 * settings_.rcrit() * rinvapprox) + (d0approx / settings_.rcrit()));
1020   int phicrit = iphi0 - 2 * irinv - 2 * id0;
1021 
1022   int iphicritmincut = settings_.phicritminmc() / globals_->ITC_L1L2()->phi0_final.K();
1023   int iphicritmaxcut = settings_.phicritmaxmc() / globals_->ITC_L1L2()->phi0_final.K();
1024 
1025   bool keepapprox = (phicritapprox > settings_.phicritminmc()) && (phicritapprox < settings_.phicritmaxmc()),
1026        keep = (phicrit > iphicritmincut) && (phicrit < iphicritmaxcut);
1027 
1028   if (settings_.debugTracklet())
1029     if (keep && !keepapprox)
1030       edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced::DDLSeeding tracklet kept with exact phicrit cut "
1031                                       "but not approximate, phicritapprox: "
1032                                    << phicritapprox;
1033   if (settings_.usephicritapprox()) {
1034     if (!keepapprox)
1035       return false;
1036   } else {
1037     if (!keep)
1038       return false;
1039   }
1040 
1041   Projection projs[N_LAYER + N_DISK];
1042 
1043   for (unsigned int i = 0; i < toR_.size(); ++i) {
1044     iphiproj[i] = phiprojapprox[i] / kphiproj;
1045     izproj[i] = zprojapprox[i] / kzproj;
1046 
1047     iphider[i] = phiderapprox[i] / kphider;
1048     izder[i] = zderapprox[i] / kzder;
1049 
1050     //check that z projection in range
1051     if (izproj[i] < -(1 << (settings_.nzbitsstub(0) - 1)))
1052       continue;
1053     if (izproj[i] >= (1 << (settings_.nzbitsstub(0) - 1)))
1054       continue;
1055 
1056     //check that phi projection in range
1057     if (iphiproj[i] >= (1 << settings_.nphibitsstub(5)) - 1)
1058       continue;
1059     if (iphiproj[i] <= 0)
1060       continue;
1061 
1062     if (rproj_[i] < settings_.rPS2S()) {
1063       iphiproj[i] >>= (settings_.nphibitsstub(5) - settings_.nphibitsstub(0));
1064     } else {
1065       izproj[i] >>= (settings_.nzbitsstub(0) - settings_.nzbitsstub(5));
1066     }
1067 
1068     if (rproj_[i] < settings_.rPS2S()) {
1069       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL123() - 1)))
1070         iphider[i] = -(1 << (settings_.nbitsphiprojderL123() - 1));
1071       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL123() - 1)))
1072         iphider[i] = (1 << (settings_.nbitsphiprojderL123() - 1)) - 1;
1073     } else {
1074       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL456() - 1)))
1075         iphider[i] = -(1 << (settings_.nbitsphiprojderL456() - 1));
1076       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL456() - 1)))
1077         iphider[i] = (1 << (settings_.nbitsphiprojderL456() - 1)) - 1;
1078     }
1079 
1080     projs[lproj_[i] - 1].init(settings_,
1081                               lproj_[i] - 1,
1082                               iphiproj[i],
1083                               izproj[i],
1084                               iphider[i],
1085                               izder[i],
1086                               phiproj[i],
1087                               zproj[i],
1088                               phider[i],
1089                               zder[i],
1090                               phiprojapprox[i],
1091                               zprojapprox[i],
1092                               phiderapprox[i],
1093                               zderapprox[i],
1094                               false);
1095   }
1096 
1097   if (std::abs(it * kt) > 1.0) {
1098     for (unsigned int i = 0; i < toZ_.size(); ++i) {
1099       iphiprojdisk[i] = phiprojdiskapprox[i] / kphiprojdisk;
1100       irprojdisk[i] = rprojdiskapprox[i] / krprojdisk;
1101 
1102       iphiderdisk[i] = phiderdiskapprox[i] / kphiderdisk;
1103       irderdisk[i] = rderdiskapprox[i] / krderdisk;
1104 
1105       if (iphiprojdisk[i] <= 0)
1106         continue;
1107       if (iphiprojdisk[i] >= (1 << settings_.nphibitsstub(0)) - 1)
1108         continue;
1109 
1110       if (irprojdisk[i] < settings_.rmindisk() / krprojdisk || irprojdisk[i] >= settings_.rmaxdisk() / krprojdisk)
1111         continue;
1112 
1113       projs[N_LAYER + i + 2].init(settings_,
1114                                   N_LAYER + i + 2,
1115                                   iphiprojdisk[i],
1116                                   irprojdisk[i],
1117                                   iphiderdisk[i],
1118                                   irderdisk[i],
1119                                   phiprojdisk[i],
1120                                   rprojdisk[i],
1121                                   phiderdisk[i],
1122                                   rderdisk[i],
1123                                   phiprojdiskapprox[i],
1124                                   rprojdiskapprox[i],
1125                                   phiderdisk[i],
1126                                   rderdisk[i],
1127                                   false);
1128     }
1129   }
1130 
1131   if (settings_.writeMonitorData("TrackletPars")) {
1132     globals_->ofstream("trackletpars.txt")
1133         << layer_ << " , " << rinv << " , " << rinvapprox << " , " << phi0 << " , " << phi0approx << " , " << t << " , "
1134         << tapprox << " , " << z0 << " , " << z0approx << " , " << d0 << " , " << d0approx << endl;
1135   }
1136 
1137   Tracklet* tracklet = new Tracklet(settings_,
1138                                     iSeed_,
1139                                     innerFPGAStub,
1140                                     middleFPGAStub,
1141                                     outerFPGAStub,
1142                                     rinv,
1143                                     phi0,
1144                                     d0,
1145                                     z0,
1146                                     t,
1147                                     rinvapprox,
1148                                     phi0approx,
1149                                     d0approx,
1150                                     z0approx,
1151                                     tapprox,
1152                                     irinv,
1153                                     iphi0,
1154                                     id0,
1155                                     iz0,
1156                                     it,
1157                                     projs,
1158                                     true);
1159 
1160   if (settings_.debugTracklet())
1161     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName()
1162                                  << " Found DDL tracklet in sector = " << iSector_ << " phi0 = " << phi0;
1163 
1164   tracklet->setTrackletIndex(trackletpars_->nTracklets());
1165   tracklet->setTCIndex(TCIndex_);
1166 
1167   if (settings_.writeMonitorData("Seeds")) {
1168     ofstream fout("seeds.txt", ofstream::app);
1169     fout << __FILE__ << ":" << __LINE__ << " " << name_ << "_" << iSector_ << " " << tracklet->getISeed() << endl;
1170     fout.close();
1171   }
1172   trackletpars_->addTracklet(tracklet);
1173 
1174   for (unsigned int j = 0; j < toR_.size(); j++) {
1175     if (settings_.debugTracklet())
1176       edm::LogVerbatim("Tracklet") << "adding layer projection " << j << "/" << toR_.size() << " " << lproj_[j] << " "
1177                                    << tracklet->validProj(lproj_[j] - 1);
1178     if (tracklet->validProj(lproj_[j] - 1)) {
1179       addLayerProj(tracklet, lproj_[j]);
1180     }
1181   }
1182 
1183   for (unsigned int j = 0; j < toZ_.size(); j++) {
1184     int disk = dproj_[j];
1185     if (disk == 0)
1186       continue;
1187     if (it < 0)
1188       disk = -disk;
1189     if (settings_.debugTracklet())
1190       edm::LogVerbatim("Tracklet") << "adding disk projection " << j << "/" << toZ_.size() << " " << disk << " "
1191                                    << tracklet->validProj(N_LAYER + abs(disk) - 1);
1192     if (tracklet->validProj(N_LAYER + abs(disk) - 1)) {
1193       addDiskProj(tracklet, disk);
1194     }
1195   }
1196 
1197   return true;
1198 }
1199 
1200 bool TrackletCalculatorDisplaced::LLDSeeding(const Stub* innerFPGAStub,
1201                                              const L1TStub* innerStub,
1202                                              const Stub* middleFPGAStub,
1203                                              const L1TStub* middleStub,
1204                                              const Stub* outerFPGAStub,
1205                                              const L1TStub* outerStub) {
1206   if (settings_.debugTracklet())
1207     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName() << " " << layer_
1208                                  << " trying stub triplet in  (L2L3D1): " << middleFPGAStub->layer().value() << " "
1209                                  << outerFPGAStub->layer().value() << " " << innerFPGAStub->disk().value();
1210 
1211   int take3 = 0;  //L2L3D1
1212   unsigned ndisks = 1;
1213 
1214   double r3 = innerStub->r();
1215   double z3 = innerStub->z();
1216   double phi3 = innerStub->phi();
1217 
1218   double r1 = middleStub->r();
1219   double z1 = middleStub->z();
1220   double phi1 = middleStub->phi();
1221 
1222   double r2 = outerStub->r();
1223   double z2 = outerStub->z();
1224   double phi2 = outerStub->phi();
1225 
1226   double rinv, phi0, d0, t, z0;
1227 
1228   double phiproj[N_LAYER - 2], zproj[N_LAYER - 2], phider[N_LAYER - 2], zder[N_LAYER - 2];
1229   double phiprojdisk[N_DISK], rprojdisk[N_DISK], phiderdisk[N_DISK], rderdisk[N_DISK];
1230 
1231   exacttracklet(r1,
1232                 z1,
1233                 phi1,
1234                 r2,
1235                 z2,
1236                 phi2,
1237                 r3,
1238                 z3,
1239                 phi3,
1240                 take3,
1241                 rinv,
1242                 phi0,
1243                 d0,
1244                 t,
1245                 z0,
1246                 phiproj,
1247                 zproj,
1248                 phiprojdisk,
1249                 rprojdisk,
1250                 phider,
1251                 zder,
1252                 phiderdisk,
1253                 rderdisk);
1254 
1255   if (settings_.debugTracklet())
1256     edm::LogVerbatim("Tracklet") << __LINE__ << ":" << __FILE__ << "LLD Exact values " << innerFPGAStub->isBarrel()
1257                                  << middleFPGAStub->isBarrel() << outerFPGAStub->isBarrel() << " " << phi3 << ", " << z3
1258                                  << ", " << r3 << ", " << phi1 << ", " << z1 << ", " << r1 << ", " << phi2 << ", " << z2
1259                                  << ", " << r2 << endl;
1260 
1261   if (settings_.useapprox()) {
1262     phi3 = innerFPGAStub->phiapprox(phimin_, phimax_);
1263     z3 = innerFPGAStub->zapprox();
1264     r3 = innerFPGAStub->rapprox();
1265 
1266     phi1 = middleFPGAStub->phiapprox(phimin_, phimax_);
1267     z1 = middleFPGAStub->zapprox();
1268     r1 = middleFPGAStub->rapprox();
1269 
1270     phi2 = outerFPGAStub->phiapprox(phimin_, phimax_);
1271     z2 = outerFPGAStub->zapprox();
1272     r2 = outerFPGAStub->rapprox();
1273   }
1274 
1275   if (settings_.debugTracklet())
1276     edm::LogVerbatim("Tracklet") << __LINE__ << ":" << __FILE__ << "LLD approx values " << innerFPGAStub->isBarrel()
1277                                  << middleFPGAStub->isBarrel() << outerFPGAStub->isBarrel() << " " << phi3 << ", " << z3
1278                                  << ", " << r3 << ", " << phi1 << ", " << z1 << ", " << r1 << ", " << phi2 << ", " << z2
1279                                  << ", " << r2 << endl;
1280 
1281   double rinvapprox, phi0approx, d0approx, tapprox, z0approx;
1282   double phiprojapprox[N_LAYER - 2], zprojapprox[N_LAYER - 2], phiderapprox[N_LAYER - 2], zderapprox[N_LAYER - 2];
1283   double phiprojdiskapprox[N_DISK], rprojdiskapprox[N_DISK];
1284   double phiderdiskapprox[N_DISK], rderdiskapprox[N_DISK];
1285 
1286   //TODO: implement the actual integer calculation
1287   if (settings_.useapprox()) {
1288     approxtracklet(r1,
1289                    z1,
1290                    phi1,
1291                    r2,
1292                    z2,
1293                    phi2,
1294                    r3,
1295                    z3,
1296                    phi3,
1297                    take3,
1298                    ndisks,
1299                    rinvapprox,
1300                    phi0approx,
1301                    d0approx,
1302                    tapprox,
1303                    z0approx,
1304                    phiprojapprox,
1305                    zprojapprox,
1306                    phiderapprox,
1307                    zderapprox,
1308                    phiprojdiskapprox,
1309                    rprojdiskapprox,
1310                    phiderdiskapprox,
1311                    rderdiskapprox);
1312   } else {
1313     rinvapprox = rinv;
1314     phi0approx = phi0;
1315     d0approx = d0;
1316     tapprox = t;
1317     z0approx = z0;
1318 
1319     for (unsigned int i = 0; i < toR_.size(); ++i) {
1320       phiprojapprox[i] = phiproj[i];
1321       zprojapprox[i] = zproj[i];
1322       phiderapprox[i] = phider[i];
1323       zderapprox[i] = zder[i];
1324     }
1325 
1326     for (unsigned int i = 0; i < toZ_.size(); ++i) {
1327       phiprojdiskapprox[i] = phiprojdisk[i];
1328       rprojdiskapprox[i] = rprojdisk[i];
1329       phiderdiskapprox[i] = phiderdisk[i];
1330       rderdiskapprox[i] = rderdisk[i];
1331     }
1332   }
1333 
1334   //store the approcximate results
1335   if (settings_.debugTracklet()) {
1336     edm::LogVerbatim("Tracklet") << "rinvapprox: " << rinvapprox << " rinv: " << rinv << endl;
1337     edm::LogVerbatim("Tracklet") << "phi0approx: " << phi0approx << " phi0: " << phi0 << endl;
1338     edm::LogVerbatim("Tracklet") << "d0approx: " << d0approx << " d0: " << d0 << endl;
1339     edm::LogVerbatim("Tracklet") << "tapprox: " << tapprox << " t: " << t << endl;
1340     edm::LogVerbatim("Tracklet") << "z0approx: " << z0approx << " z0: " << z0 << endl;
1341   }
1342 
1343   for (unsigned int i = 0; i < toR_.size(); ++i) {
1344     if (settings_.debugTracklet()) {
1345       edm::LogVerbatim("Tracklet") << "phiprojapprox[" << i << "]: " << phiprojapprox[i] << " phiproj[" << i
1346                                    << "]: " << phiproj[i] << endl;
1347       edm::LogVerbatim("Tracklet") << "zprojapprox[" << i << "]: " << zprojapprox[i] << " zproj[" << i
1348                                    << "]: " << zproj[i] << endl;
1349       edm::LogVerbatim("Tracklet") << "phiderapprox[" << i << "]: " << phiderapprox[i] << " phider[" << i
1350                                    << "]: " << phider[i] << endl;
1351       edm::LogVerbatim("Tracklet") << "zderapprox[" << i << "]: " << zderapprox[i] << " zder[" << i << "]: " << zder[i]
1352                                    << endl;
1353     }
1354   }
1355 
1356   for (unsigned int i = 0; i < toZ_.size(); ++i) {
1357     if (settings_.debugTracklet()) {
1358       edm::LogVerbatim("Tracklet") << "phiprojdiskapprox[" << i << "]: " << phiprojdiskapprox[i] << " phiprojdisk[" << i
1359                                    << "]: " << phiprojdisk[i] << endl;
1360       edm::LogVerbatim("Tracklet") << "rprojdiskapprox[" << i << "]: " << rprojdiskapprox[i] << " rprojdisk[" << i
1361                                    << "]: " << rprojdisk[i] << endl;
1362       edm::LogVerbatim("Tracklet") << "phiderdiskapprox[" << i << "]: " << phiderdiskapprox[i] << " phiderdisk[" << i
1363                                    << "]: " << phiderdisk[i] << endl;
1364       edm::LogVerbatim("Tracklet") << "rderdiskapprox[" << i << "]: " << rderdiskapprox[i] << " rderdisk[" << i
1365                                    << "]: " << rderdisk[i] << endl;
1366     }
1367   }
1368 
1369   //now binary
1370   double krinv = settings_.kphi1() / settings_.kr() * pow(2, settings_.rinv_shift()),
1371          kphi0 = settings_.kphi1() * pow(2, settings_.phi0_shift()),
1372          kt = settings_.kz() / settings_.kr() * pow(2, settings_.t_shift()),
1373          kz0 = settings_.kz() * pow(2, settings_.z0_shift()),
1374          kphiproj = settings_.kphi1() * pow(2, settings_.SS_phiL_shift()),
1375          kphider = settings_.kphi1() / settings_.kr() * pow(2, settings_.SS_phiderL_shift()),
1376          kzproj = settings_.kz() * pow(2, settings_.PS_zL_shift()),
1377          kzder = settings_.kz() / settings_.kr() * pow(2, settings_.PS_zderL_shift()),
1378          kphiprojdisk = settings_.kphi1() * pow(2, settings_.SS_phiD_shift()),
1379          kphiderdisk = settings_.kphi1() / settings_.kr() * pow(2, settings_.SS_phiderD_shift()),
1380          krprojdisk = settings_.kr() * pow(2, settings_.PS_rD_shift()),
1381          krderdisk = settings_.kr() / settings_.kz() * pow(2, settings_.PS_rderD_shift());
1382 
1383   int irinv, iphi0, id0, it, iz0;
1384   int iphiproj[N_LAYER - 2], izproj[N_LAYER - 2], iphider[N_LAYER - 2], izder[N_LAYER - 2];
1385   int iphiprojdisk[N_DISK], irprojdisk[N_DISK], iphiderdisk[N_DISK], irderdisk[N_DISK];
1386 
1387   //store the binary results
1388   irinv = rinvapprox / krinv;
1389   iphi0 = phi0approx / kphi0;
1390   id0 = d0approx / settings_.kd0();
1391   it = tapprox / kt;
1392   iz0 = z0approx / kz0;
1393 
1394   bool success = true;
1395   if (std::abs(rinvapprox) > settings_.rinvcut()) {
1396     if (settings_.debugTracklet())
1397       edm::LogVerbatim("Tracklet") << "TrackletCalculator:: LLD Seeding irinv too large: " << rinvapprox << "(" << irinv
1398                                    << ")";
1399     success = false;
1400   }
1401   if (std::abs(z0approx) > settings_.disp_z0cut()) {
1402     if (settings_.debugTracklet())
1403       edm::LogVerbatim("Tracklet") << "Failed tracklet z0 cut " << z0approx;
1404     success = false;
1405   }
1406   if (std::abs(d0approx) > settings_.maxd0()) {
1407     if (settings_.debugTracklet())
1408       edm::LogVerbatim("Tracklet") << "Failed tracklet approx d0 cut " << d0approx;
1409     success = false;
1410   }
1411   if (std::abs(d0) > settings_.maxd0()) {
1412     if (settings_.debugTracklet())
1413       edm::LogVerbatim("Tracklet") << "Failed tracklet exact d0 cut " << d0;
1414     success = false;
1415   }
1416 
1417   if (!success)
1418     return false;
1419 
1420   double phicritapprox = phi0approx - asin((0.5 * settings_.rcrit() * rinvapprox) + (d0approx / settings_.rcrit()));
1421   int phicrit = iphi0 - 2 * irinv - 2 * id0;
1422 
1423   int iphicritmincut = settings_.phicritminmc() / globals_->ITC_L1L2()->phi0_final.K();
1424   int iphicritmaxcut = settings_.phicritmaxmc() / globals_->ITC_L1L2()->phi0_final.K();
1425 
1426   bool keepapprox = (phicritapprox > settings_.phicritminmc()) && (phicritapprox < settings_.phicritmaxmc()),
1427        keep = (phicrit > iphicritmincut) && (phicrit < iphicritmaxcut);
1428 
1429   if (settings_.debugTracklet())
1430     if (keep && !keepapprox)
1431       edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced::LLDSeeding tracklet kept with exact phicrit cut "
1432                                       "but not approximate, phicritapprox: "
1433                                    << phicritapprox;
1434   if (settings_.usephicritapprox()) {
1435     if (!keepapprox)
1436       return false;
1437   } else {
1438     if (!keep)
1439       return false;
1440   }
1441 
1442   Projection projs[N_LAYER + N_DISK];
1443 
1444   for (unsigned int i = 0; i < toR_.size(); ++i) {
1445     iphiproj[i] = phiprojapprox[i] / kphiproj;
1446     izproj[i] = zprojapprox[i] / kzproj;
1447 
1448     iphider[i] = phiderapprox[i] / kphider;
1449     izder[i] = zderapprox[i] / kzder;
1450 
1451     if (izproj[i] < -(1 << (settings_.nzbitsstub(0) - 1)))
1452       continue;
1453     if (izproj[i] >= (1 << (settings_.nzbitsstub(0) - 1)))
1454       continue;
1455 
1456     //this is left from the original....
1457     if (iphiproj[i] >= (1 << settings_.nphibitsstub(5)) - 1)
1458       continue;
1459     if (iphiproj[i] <= 0)
1460       continue;
1461 
1462     if (rproj_[i] < settings_.rPS2S()) {
1463       iphiproj[i] >>= (settings_.nphibitsstub(5) - settings_.nphibitsstub(0));
1464     } else {
1465       izproj[i] >>= (settings_.nzbitsstub(0) - settings_.nzbitsstub(5));
1466     }
1467 
1468     if (rproj_[i] < settings_.rPS2S()) {
1469       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL123() - 1)))
1470         iphider[i] = -(1 << (settings_.nbitsphiprojderL123() - 1));
1471       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL123() - 1)))
1472         iphider[i] = (1 << (settings_.nbitsphiprojderL123() - 1)) - 1;
1473     } else {
1474       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL456() - 1)))
1475         iphider[i] = -(1 << (settings_.nbitsphiprojderL456() - 1));
1476       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL456() - 1)))
1477         iphider[i] = (1 << (settings_.nbitsphiprojderL456() - 1)) - 1;
1478     }
1479 
1480     projs[lproj_[i] - 1].init(settings_,
1481                               lproj_[i] - 1,
1482                               iphiproj[i],
1483                               izproj[i],
1484                               iphider[i],
1485                               izder[i],
1486                               phiproj[i],
1487                               zproj[i],
1488                               phider[i],
1489                               zder[i],
1490                               phiprojapprox[i],
1491                               zprojapprox[i],
1492                               phiderapprox[i],
1493                               zderapprox[i],
1494                               false);
1495   }
1496 
1497   if (std::abs(it * kt) > 1.0) {
1498     for (unsigned int i = 0; i < toZ_.size(); ++i) {
1499       iphiprojdisk[i] = phiprojdiskapprox[i] / kphiprojdisk;
1500       irprojdisk[i] = rprojdiskapprox[i] / krprojdisk;
1501 
1502       iphiderdisk[i] = phiderdiskapprox[i] / kphiderdisk;
1503       irderdisk[i] = rderdiskapprox[i] / krderdisk;
1504 
1505       //Check phi range of projection
1506       if (iphiprojdisk[i] <= 0)
1507         continue;
1508       if (iphiprojdisk[i] >= (1 << settings_.nphibitsstub(0)) - 1)
1509         continue;
1510 
1511       //Check r range of projection
1512       if (irprojdisk[i] < settings_.rmindisk() / krprojdisk || irprojdisk[i] >= settings_.rmaxdisk() / krprojdisk)
1513         continue;
1514 
1515       projs[N_LAYER + i + 1].init(settings_,
1516                                   N_LAYER + i + 1,
1517                                   iphiprojdisk[i],
1518                                   irprojdisk[i],
1519                                   iphiderdisk[i],
1520                                   irderdisk[i],
1521                                   phiprojdisk[i],
1522                                   rprojdisk[i],
1523                                   phiderdisk[i],
1524                                   rderdisk[i],
1525                                   phiprojdiskapprox[i],
1526                                   rprojdiskapprox[i],
1527                                   phiderdisk[i],
1528                                   rderdisk[i],
1529                                   false);
1530     }
1531   }
1532 
1533   if (settings_.writeMonitorData("TrackletPars")) {
1534     globals_->ofstream("trackletpars.txt")
1535         << layer_ << " , " << rinv << " , " << rinvapprox << " , " << phi0 << " , " << phi0approx << " , " << t << " , "
1536         << tapprox << " , " << z0 << " , " << z0approx << " , " << d0 << " , " << d0approx << endl;
1537   }
1538 
1539   Tracklet* tracklet = new Tracklet(settings_,
1540                                     iSeed_,
1541                                     innerFPGAStub,
1542                                     middleFPGAStub,
1543                                     outerFPGAStub,
1544                                     rinv,
1545                                     phi0,
1546                                     d0,
1547                                     z0,
1548                                     t,
1549                                     rinvapprox,
1550                                     phi0approx,
1551                                     d0approx,
1552                                     z0approx,
1553                                     tapprox,
1554                                     irinv,
1555                                     iphi0,
1556                                     id0,
1557                                     iz0,
1558                                     it,
1559                                     projs,
1560                                     false);
1561 
1562   if (settings_.debugTracklet())
1563     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName()
1564                                  << " Found LLD tracklet in sector = " << iSector_ << " phi0 = " << phi0;
1565 
1566   tracklet->setTrackletIndex(trackletpars_->nTracklets());
1567   tracklet->setTCIndex(TCIndex_);
1568 
1569   if (settings_.writeMonitorData("Seeds")) {
1570     ofstream fout("seeds.txt", ofstream::app);
1571     fout << __FILE__ << ":" << __LINE__ << " " << name_ << "_" << iSector_ << " " << tracklet->getISeed() << endl;
1572     fout.close();
1573   }
1574   trackletpars_->addTracklet(tracklet);
1575 
1576   for (unsigned int j = 0; j < toR_.size(); j++) {
1577     if (settings_.debugTracklet())
1578       edm::LogVerbatim("Tracklet") << "adding layer projection " << j << "/" << toR_.size() << " " << lproj_[j];
1579     if (tracklet->validProj(lproj_[j] - 1)) {
1580       addLayerProj(tracklet, lproj_[j]);
1581     }
1582   }
1583 
1584   for (unsigned int j = 0; j < toZ_.size(); j++) {
1585     int disk = dproj_[j];
1586     if (disk == 0)
1587       continue;
1588     if (it < 0)
1589       disk = -disk;
1590     if (settings_.debugTracklet())
1591       edm::LogVerbatim("Tracklet") << "adding disk projection " << j << "/" << toZ_.size() << " " << disk;
1592     if (tracklet->validProj(N_LAYER + abs(disk) - 1)) {
1593       addDiskProj(tracklet, disk);
1594     }
1595   }
1596 
1597   return true;
1598 }
1599 
1600 void TrackletCalculatorDisplaced::exactproj(double rproj,
1601                                             double rinv,
1602                                             double phi0,
1603                                             double d0,
1604                                             double t,
1605                                             double z0,
1606                                             double r0,
1607                                             double& phiproj,
1608                                             double& zproj,
1609                                             double& phider,
1610                                             double& zder) {
1611   double rho = 1 / rinv;
1612   if (rho < 0) {
1613     r0 = -r0;
1614   }
1615   phiproj = phi0 - asin((rproj * rproj + r0 * r0 - rho * rho) / (2 * rproj * r0));
1616   double beta = acos((rho * rho + r0 * r0 - rproj * rproj) / (2 * r0 * rho));
1617   zproj = z0 + t * std::abs(rho * beta);
1618 
1619   //not exact, but close
1620   phider = -0.5 * rinv / sqrt(1 - pow(0.5 * rproj * rinv, 2)) + d0 / (rproj * rproj);
1621   zder = t / sqrt(1 - pow(0.5 * rproj * rinv, 2));
1622 
1623   if (settings_.debugTracklet())
1624     edm::LogVerbatim("Tracklet") << "exact proj layer at " << rproj << " : " << phiproj << " " << zproj;
1625 }
1626 
1627 void TrackletCalculatorDisplaced::exactprojdisk(double zproj,
1628                                                 double rinv,
1629                                                 double,
1630                                                 double,  //phi0 and d0 are not used.
1631                                                 double t,
1632                                                 double z0,
1633                                                 double x0,
1634                                                 double y0,
1635                                                 double& phiproj,
1636                                                 double& rproj,
1637                                                 double& phider,
1638                                                 double& rder) {
1639   //protect against t=0
1640   if (std::abs(t) < 0.1)
1641     t = 0.1;
1642   if (t < 0)
1643     zproj = -zproj;
1644   double rho = std::abs(1 / rinv);
1645   double beta = (zproj - z0) / (t * rho);
1646   double phiV = atan2(-y0, -x0);
1647   double c = rinv > 0 ? -1 : 1;
1648 
1649   double x = x0 + rho * cos(phiV + c * beta);
1650   double y = y0 + rho * sin(phiV + c * beta);
1651 
1652   phiproj = atan2(y, x);
1653 
1654   phiproj = reco::reduceRange(phiproj - phimin_);
1655 
1656   rproj = sqrt(x * x + y * y);
1657 
1658   phider = c / t / (x * x + y * y) * (rho + x0 * cos(phiV + c * beta) + y0 * sin(phiV + c * beta));
1659   rder = c / t / rproj * (y0 * cos(phiV + c * beta) - x0 * sin(phiV + c * beta));
1660 
1661   if (settings_.debugTracklet())
1662     edm::LogVerbatim("Tracklet") << "exact proj disk at" << zproj << " : " << phiproj << " " << rproj;
1663 }
1664 
1665 void TrackletCalculatorDisplaced::exacttracklet(double r1,
1666                                                 double z1,
1667                                                 double phi1,
1668                                                 double r2,
1669                                                 double z2,
1670                                                 double phi2,
1671                                                 double r3,
1672                                                 double z3,
1673                                                 double phi3,
1674                                                 int take3,
1675                                                 double& rinv,
1676                                                 double& phi0,
1677                                                 double& d0,
1678                                                 double& t,
1679                                                 double& z0,
1680                                                 double phiproj[N_LAYER - 2],
1681                                                 double zproj[N_LAYER - 2],
1682                                                 double phiprojdisk[N_DISK],
1683                                                 double rprojdisk[N_DISK],
1684                                                 double phider[N_LAYER - 2],
1685                                                 double zder[N_LAYER - 2],
1686                                                 double phiderdisk[N_DISK],
1687                                                 double rderdisk[N_DISK]) {
1688   //two lines perpendicular to the 1->2 and 2->3
1689   double x1 = r1 * cos(phi1);
1690   double x2 = r2 * cos(phi2);
1691   double x3 = r3 * cos(phi3);
1692 
1693   double y1 = r1 * sin(phi1);
1694   double y2 = r2 * sin(phi2);
1695   double y3 = r3 * sin(phi3);
1696 
1697   double dy21 = y2 - y1;
1698   double dy32 = y3 - y2;
1699 
1700   //Hack to protect against dividing by zero
1701   //code should be rewritten to avoid this
1702   if (dy21 == 0.0)
1703     dy21 = 1e-9;
1704   if (dy32 == 0.0)
1705     dy32 = 1e-9;
1706 
1707   double k1 = -(x2 - x1) / dy21;
1708   double k2 = -(x3 - x2) / dy32;
1709   double b1 = 0.5 * (y2 + y1) - 0.5 * (x1 + x2) * k1;
1710   double b2 = 0.5 * (y3 + y2) - 0.5 * (x2 + x3) * k2;
1711   //their intersection gives the center of the circle
1712   double y0 = (b1 * k2 - b2 * k1) / (k2 - k1);
1713   double x0 = (b1 - b2) / (k2 - k1);
1714   //get the radius three ways:
1715   double R1 = sqrt(pow(x1 - x0, 2) + pow(y1 - y0, 2));
1716   double R2 = sqrt(pow(x2 - x0, 2) + pow(y2 - y0, 2));
1717   double R3 = sqrt(pow(x3 - x0, 2) + pow(y3 - y0, 2));
1718   //check if the same
1719   double eps1 = std::abs(R1 / R2 - 1);
1720   double eps2 = std::abs(R3 / R2 - 1);
1721   if (eps1 > 1e-10 || eps2 > 1e-10)
1722     edm::LogVerbatim("Tracklet") << "&&&&&&&&&&&& bad circle! " << R1 << "\t" << R2 << "\t" << R3;
1723 
1724   if (settings_.debugTracklet())
1725     edm::LogVerbatim("Tracklet") << "phimin_: " << phimin_ << " phimax_: " << phimax_;
1726   //results
1727   rinv = 1. / R1;
1728   phi0 = 0.5 * M_PI + atan2(y0, x0);
1729 
1730   phi0 -= phimin_;
1731 
1732   d0 = -R1 + sqrt(x0 * x0 + y0 * y0);
1733   //sign of rinv:
1734   double dphi = reco::reduceRange(phi3 - atan2(y0, x0));
1735   if (dphi < 0) {
1736     rinv = -rinv;
1737     d0 = -d0;
1738     phi0 = phi0 + M_PI;
1739   }
1740   phi0 = angle0to2pi::make0To2pi(phi0);
1741 
1742   //now in RZ:
1743   //turning angle
1744   double beta1 = reco::reduceRange(atan2(y1 - y0, x1 - x0) - atan2(-y0, -x0));
1745   double beta2 = reco::reduceRange(atan2(y2 - y0, x2 - x0) - atan2(-y0, -x0));
1746   double beta3 = reco::reduceRange(atan2(y3 - y0, x3 - x0) - atan2(-y0, -x0));
1747 
1748   double t12 = (z2 - z1) / std::abs(beta2 - beta1) / R1;
1749   double z12 = (z1 * beta2 - z2 * beta1) / (beta2 - beta1);
1750   double t13 = (z3 - z1) / std::abs(beta3 - beta1) / R1;
1751   double z13 = (z1 * beta3 - z3 * beta1) / (beta3 - beta1);
1752 
1753   if (take3 > 0) {
1754     //take 13 (large lever arm)
1755     t = t13;
1756     z0 = z13;
1757   } else {
1758     //take 12 (pixel layers)
1759     t = t12;
1760     z0 = z12;
1761   }
1762 
1763   for (unsigned int i = 0; i < toR_.size(); i++) {
1764     exactproj(toR_[i], rinv, phi0, d0, t, z0, sqrt(x0 * x0 + y0 * y0), phiproj[i], zproj[i], phider[i], zder[i]);
1765   }
1766 
1767   for (unsigned int i = 0; i < toZ_.size(); i++) {
1768     exactprojdisk(toZ_[i], rinv, phi0, d0, t, z0, x0, y0, phiprojdisk[i], rprojdisk[i], phiderdisk[i], rderdisk[i]);
1769   }
1770 
1771   if (settings_.debugTracklet())
1772     edm::LogVerbatim("Tracklet") << "exact tracklet: " << rinv << " " << phi0 << " " << t << " " << z0 << " " << d0;
1773 }
1774 
1775 void TrackletCalculatorDisplaced::approxproj(double halfRinv,
1776                                              double phi0,
1777                                              double d0,
1778                                              double t,
1779                                              double z0,
1780                                              double halfRinv_0,
1781                                              double d0_0,  // zeroeth order result for higher order terms calculation
1782                                              double rmean,
1783                                              double& phiproj,
1784                                              double& phiprojder,
1785                                              double& zproj,
1786                                              double& zprojder) {
1787   if (std::abs(2.0 * halfRinv) > settings_.rinvcut() || std::abs(z0) > settings_.disp_z0cut() ||
1788       std::abs(d0) > settings_.maxd0()) {
1789     phiproj = 0.0;
1790     return;
1791   }
1792   double rmeanInv = 1.0 / rmean;
1793 
1794   phiproj = phi0 + rmean * (-halfRinv + 2.0 * d0_0 * halfRinv_0 * halfRinv_0) +
1795             rmeanInv * (-d0 + halfRinv_0 * d0_0 * d0_0) + sixth * pow(-rmean * halfRinv_0 - rmeanInv * d0_0, 3);
1796   phiprojder = -halfRinv + d0 * rmeanInv * rmeanInv;  //removed all high terms
1797 
1798   zproj = z0 + t * rmean - 0.5 * rmeanInv * t * d0_0 * d0_0 - t * rmean * halfRinv * d0 +
1799           sixth * pow(rmean, 3) * t * halfRinv_0 * halfRinv_0;
1800   zprojder = t;  // removed all high terms
1801 
1802   phiproj = angle0to2pi::make0To2pi(phiproj);
1803 
1804   if (settings_.debugTracklet())
1805     edm::LogVerbatim("Tracklet") << "approx proj layer at " << rmean << " : " << phiproj << " " << zproj << endl;
1806 }
1807 
1808 void TrackletCalculatorDisplaced::approxprojdisk(double halfRinv,
1809                                                  double phi0,
1810                                                  double d0,
1811                                                  double t,
1812                                                  double z0,
1813                                                  double halfRinv_0,
1814                                                  double d0_0,  // zeroeth order result for higher order terms calculation
1815                                                  double zmean,
1816                                                  double& phiproj,
1817                                                  double& phiprojder,
1818                                                  double& rproj,
1819                                                  double& rprojder) {
1820   if (std::abs(2.0 * halfRinv) > settings_.rinvcut() || std::abs(z0) > settings_.disp_z0cut() ||
1821       std::abs(d0) > settings_.maxd0()) {
1822     phiproj = 0.0;
1823     return;
1824   }
1825 
1826   if (t < 0)
1827     zmean = -zmean;
1828 
1829   double zmeanInv = 1.0 / zmean, rstar = (zmean - z0) / t,
1830          epsilon = 0.5 * zmeanInv * zmeanInv * d0_0 * d0_0 * t * t + halfRinv * d0 -
1831                    sixth * rstar * rstar * halfRinv_0 * halfRinv_0;
1832 
1833   rproj = rstar * (1 + epsilon);
1834   rprojder = 1 / t;
1835 
1836   double A = rproj * halfRinv;
1837   double B = -d0 * t * zmeanInv * (1 + z0 * zmeanInv) * (1 - epsilon);
1838   double C = -d0 * halfRinv;
1839   double A_0 = rproj * halfRinv_0;
1840   double B_0 = -d0_0 * t * zmeanInv * (1 + z0 * zmeanInv) * (1 - epsilon);
1841   // double C_0 = -d0_0 * halfRinv_0;
1842 
1843   phiproj = phi0 - A + B * (1 + C - 2 * A_0 * A_0) + sixth * pow(-A_0 + B_0, 3);
1844   phiprojder = -halfRinv / t + d0 * t * zmeanInv * zmeanInv;
1845 
1846   phiproj = angle0to2pi::make0To2pi(phiproj);
1847 
1848   if (settings_.debugTracklet())
1849     edm::LogVerbatim("Tracklet") << "approx proj disk at" << zmean << " : " << phiproj << " " << rproj << endl;
1850 }
1851 
1852 void TrackletCalculatorDisplaced::approxtracklet(double r1,
1853                                                  double z1,
1854                                                  double phi1,
1855                                                  double r2,
1856                                                  double z2,
1857                                                  double phi2,
1858                                                  double r3,
1859                                                  double z3,
1860                                                  double phi3,
1861                                                  bool take3,
1862                                                  unsigned ndisks,
1863                                                  double& rinv,
1864                                                  double& phi0,
1865                                                  double& d0,
1866                                                  double& t,
1867                                                  double& z0,
1868                                                  double phiproj[4],
1869                                                  double zproj[4],
1870                                                  double phider[4],
1871                                                  double zder[4],
1872                                                  double phiprojdisk[5],
1873                                                  double rprojdisk[5],
1874                                                  double phiderdisk[5],
1875                                                  double rderdisk[5]) {
1876   double a = 1.0 / ((r1 - r2) * (r1 - r3));
1877   double b = 1.0 / ((r1 - r2) * (r2 - r3));
1878   double c = 1.0 / ((r1 - r3) * (r2 - r3));
1879 
1880   // first iteration in r-phi plane
1881   double halfRinv_0 = -phi1 * r1 * a + phi2 * r2 * b - phi3 * r3 * c;
1882   double d0_0 = r1 * r2 * r3 * (-phi1 * a + phi2 * b - phi3 * c);
1883 
1884   // corrections to phi1, phi2, and phi3
1885   double r = r2, z = z2;
1886   if (take3)
1887     r = r3, z = z3;
1888 
1889   double d0OverR1 = d0_0 * rzmeanInv_[0] * (ndisks > 2 ? std::abs((z - z1) / (r - r1)) : 1.0);
1890   double d0OverR2 = d0_0 * rzmeanInv_[1] * (ndisks > 1 ? std::abs((z - z1) / (r - r1)) : 1.0);
1891   double d0OverR3 = d0_0 * rzmeanInv_[2] * (ndisks > 0 ? std::abs((z - z1) / (r - r1)) : 1.0);
1892 
1893   double d0OverR = d0OverR2;
1894   if (take3)
1895     d0OverR = d0OverR3;
1896 
1897   double c1 = d0_0 * halfRinv_0 * d0OverR1 + 2.0 * d0_0 * halfRinv_0 * r1 * halfRinv_0 +
1898               sixth * pow(-r1 * halfRinv_0 - d0OverR1, 3);
1899   double c2 = d0_0 * halfRinv_0 * d0OverR2 + 2.0 * d0_0 * halfRinv_0 * r2 * halfRinv_0 +
1900               sixth * pow(-r2 * halfRinv_0 - d0OverR2, 3);
1901   double c3 = d0_0 * halfRinv_0 * d0OverR3 + 2.0 * d0_0 * halfRinv_0 * r3 * halfRinv_0 +
1902               sixth * pow(-r3 * halfRinv_0 - d0OverR3, 3);
1903 
1904   double phi1c = phi1 - c1;
1905   double phi2c = phi2 - c2;
1906   double phi3c = phi3 - c3;
1907 
1908   // second iteration in r-phi plane
1909   double halfRinv = -phi1c * r1 * a + phi2c * r2 * b - phi3c * r3 * c;
1910   phi0 = -phi1c * r1 * (r2 + r3) * a + phi2c * r2 * (r1 + r3) * b - phi3c * r3 * (r1 + r2) * c;
1911   d0 = r1 * r2 * r3 * (-phi1c * a + phi2c * b - phi3c * c);
1912 
1913   t = ((z - z1) / (r - r1)) *
1914       (1. + d0 * halfRinv - 0.5 * d0OverR1 * d0OverR - sixth * (r1 * r1 + r2 * r2 + r1 * r2) * halfRinv_0 * halfRinv_0);
1915   z0 = z1 - t * r1 * (1.0 - d0_0 * halfRinv_0 - 0.5 * d0OverR1 * d0OverR1 + sixth * r1 * r1 * halfRinv_0 * halfRinv_0);
1916 
1917   rinv = 2.0 * halfRinv;
1918   phi0 += -phimin_;
1919 
1920   phi0 = angle0to2pi::make0To2pi(phi0);
1921 
1922   for (unsigned int i = 0; i < toR_.size(); i++) {
1923     approxproj(halfRinv,
1924                phi0,
1925                d0,
1926                t,
1927                z0,
1928                halfRinv_0,
1929                d0_0,  // added _0 version for high term calculations
1930                toR_.at(i),
1931                phiproj[i],
1932                phider[i],
1933                zproj[i],
1934                zder[i]);
1935   }
1936 
1937   for (unsigned int i = 0; i < toZ_.size(); i++) {
1938     approxprojdisk(halfRinv,
1939                    phi0,
1940                    d0,
1941                    t,
1942                    z0,
1943                    halfRinv_0,
1944                    d0_0,  // added _0 version for high term calculations
1945                    toZ_.at(i),
1946                    phiprojdisk[i],
1947                    phiderdisk[i],
1948                    rprojdisk[i],
1949                    rderdisk[i]);
1950   }
1951 
1952   if (std::abs(rinv) > settings_.rinvcut() || std::abs(z0) > settings_.disp_z0cut() ||
1953       std::abs(d0) > settings_.maxd0()) {
1954     phi0 = 0.0;
1955     return;
1956   }
1957 
1958   if (settings_.debugTracklet())
1959     edm::LogVerbatim("Tracklet") << "TCD approx tracklet: " << rinv << " " << phi0 << " " << t << " " << z0 << " " << d0
1960                                  << endl;
1961 }