Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2021-07-07 22:33:38 

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 d0 cut " << d0approx;
0587     success = false;
0588   }
0589 
0590   if (!success) {
0591     return false;
0592   }
0593 
0594   double phicritapprox = phi0approx - asin(0.5 * settings_.rcrit() * rinvapprox);
0595   int phicrit = iphi0 - 2 * irinv;
0596 
0597   int iphicritmincut = settings_.phicritminmc() / globals_->ITC_L1L2()->phi0_final.K();
0598   int iphicritmaxcut = settings_.phicritmaxmc() / globals_->ITC_L1L2()->phi0_final.K();
0599 
0600   bool keepapprox = (phicritapprox > settings_.phicritminmc()) && (phicritapprox < settings_.phicritmaxmc()),
0601        keep = (phicrit > iphicritmincut) && (phicrit < iphicritmaxcut);
0602 
0603   if (settings_.debugTracklet())
0604     if (keep && !keepapprox)
0605       edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced::LLLSeeding tracklet kept with exact phicrit cut "
0606                                       "but not approximate, phicritapprox: "
0607                                    << phicritapprox;
0608   if (settings_.usephicritapprox()) {
0609     if (!keepapprox) {
0610       return false;
0611     }
0612   } else {
0613     if (!keep) {
0614       return false;
0615     }
0616   }
0617 
0618   for (unsigned int i = 0; i < toR_.size(); ++i) {
0619     iphiproj[i] = phiprojapprox[i] / kphiproj;
0620     izproj[i] = zprojapprox[i] / kzproj;
0621 
0622     iphider[i] = phiderapprox[i] / kphider;
0623     izder[i] = zderapprox[i] / kzder;
0624 
0625     //check that z projection is in range
0626     if (izproj[i] < -(1 << (settings_.nzbitsstub(0) - 1)))
0627       continue;
0628     if (izproj[i] >= (1 << (settings_.nzbitsstub(0) - 1)))
0629       continue;
0630 
0631     //check that phi projection is in range
0632     if (iphiproj[i] >= (1 << settings_.nphibitsstub(N_LAYER - 1)) - 1)
0633       continue;
0634     if (iphiproj[i] <= 0)
0635       continue;
0636 
0637     //adjust number of bits for phi and z projection
0638     if (rproj_[i] < settings_.rPS2S()) {
0639       iphiproj[i] >>= (settings_.nphibitsstub(N_LAYER - 1) - settings_.nphibitsstub(0));
0640       if (iphiproj[i] >= (1 << settings_.nphibitsstub(0)) - 1)
0641         iphiproj[i] = (1 << settings_.nphibitsstub(0)) - 2;  //-2 not to hit atExtreme
0642     } else {
0643       izproj[i] >>= (settings_.nzbitsstub(0) - settings_.nzbitsstub(N_LAYER - 1));
0644     }
0645 
0646     if (rproj_[i] < settings_.rPS2S()) {
0647       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL123() - 1))) {
0648         iphider[i] = -(1 << (settings_.nbitsphiprojderL123() - 1));
0649       }
0650       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL123() - 1))) {
0651         iphider[i] = (1 << (settings_.nbitsphiprojderL123() - 1)) - 1;
0652       }
0653     } else {
0654       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL456() - 1))) {
0655         iphider[i] = -(1 << (settings_.nbitsphiprojderL456() - 1));
0656       }
0657       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL456() - 1))) {
0658         iphider[i] = (1 << (settings_.nbitsphiprojderL456() - 1)) - 1;
0659       }
0660     }
0661 
0662     projs[lproj_[i] - 1].init(settings_,
0663                               lproj_[i] - 1,
0664                               iphiproj[i],
0665                               izproj[i],
0666                               iphider[i],
0667                               izder[i],
0668                               phiproj[i],
0669                               zproj[i],
0670                               phider[i],
0671                               zder[i],
0672                               phiprojapprox[i],
0673                               zprojapprox[i],
0674                               phiderapprox[i],
0675                               zderapprox[i],
0676                               false);
0677   }
0678 
0679   if (std::abs(it * kt) > 1.0) {
0680     for (unsigned int i = 0; i < toZ_.size(); ++i) {
0681       iphiprojdisk[i] = phiprojdiskapprox[i] / kphiprojdisk;
0682       irprojdisk[i] = rprojdiskapprox[i] / krprojdisk;
0683 
0684       iphiderdisk[i] = phiderdiskapprox[i] / kphiderdisk;
0685       irderdisk[i] = rderdiskapprox[i] / krderdisk;
0686 
0687       //check phi projection in range
0688       if (iphiprojdisk[i] <= 0)
0689         continue;
0690       if (iphiprojdisk[i] >= (1 << settings_.nphibitsstub(0)) - 1)
0691         continue;
0692 
0693       //check r projection in range
0694       if (rprojdiskapprox[i] < settings_.rmindisk() || rprojdiskapprox[i] > settings_.rmaxdisk())
0695         continue;
0696 
0697       projs[N_LAYER + i].init(settings_,
0698                               N_LAYER + i,
0699                               iphiprojdisk[i],
0700                               irprojdisk[i],
0701                               iphiderdisk[i],
0702                               irderdisk[i],
0703                               phiprojdisk[i],
0704                               rprojdisk[i],
0705                               phiderdisk[i],
0706                               rderdisk[i],
0707                               phiprojdiskapprox[i],
0708                               rprojdiskapprox[i],
0709                               phiderdisk[i],
0710                               rderdisk[i],
0711                               false);
0712     }
0713   }
0714 
0715   if (settings_.writeMonitorData("TrackletPars")) {
0716     globals_->ofstream("trackletpars.txt")
0717         << layer_ << " , " << rinv << " , " << rinvapprox << " , " << phi0 << " , " << phi0approx << " , " << t << " , "
0718         << tapprox << " , " << z0 << " , " << z0approx << " , " << d0 << " , " << d0approx << endl;
0719   }
0720 
0721   Tracklet* tracklet = new Tracklet(settings_,
0722                                     iSeed_,
0723                                     innerFPGAStub,
0724                                     middleFPGAStub,
0725                                     outerFPGAStub,
0726                                     rinv,
0727                                     phi0,
0728                                     d0,
0729                                     z0,
0730                                     t,
0731                                     rinvapprox,
0732                                     phi0approx,
0733                                     d0approx,
0734                                     z0approx,
0735                                     tapprox,
0736                                     irinv,
0737                                     iphi0,
0738                                     id0,
0739                                     iz0,
0740                                     it,
0741                                     projs,
0742                                     false);
0743 
0744   if (settings_.debugTracklet())
0745     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName()
0746                                  << " Found LLL tracklet in sector = " << iSector_ << " phi0 = " << phi0;
0747 
0748   tracklet->setTrackletIndex(trackletpars_->nTracklets());
0749   tracklet->setTCIndex(TCIndex_);
0750 
0751   if (settings_.writeMonitorData("Seeds")) {
0752     ofstream fout("seeds.txt", ofstream::app);
0753     fout << __FILE__ << ":" << __LINE__ << " " << name_ << "_" << iSector_ << " " << tracklet->getISeed() << endl;
0754     fout.close();
0755   }
0756   trackletpars_->addTracklet(tracklet);
0757 
0758   bool addL5 = false;
0759   bool addL6 = false;
0760   for (unsigned int j = 0; j < toR_.size(); j++) {
0761     bool added = false;
0762 
0763     if (settings_.debugTracklet())
0764       edm::LogVerbatim("Tracklet") << "adding layer projection " << j << "/" << toR_.size() << " " << lproj_[j];
0765     if (tracklet->validProj(lproj_[j] - 1)) {
0766       added = addLayerProj(tracklet, lproj_[j]);
0767       if (added && lproj_[j] == 5)
0768         addL5 = true;
0769       if (added && lproj_[j] == 6)
0770         addL6 = true;
0771     }
0772   }
0773 
0774   for (unsigned int j = 0; j < toZ_.size(); j++) {
0775     int disk = dproj_[j];
0776     if (disk == 0)
0777       continue;
0778     if (disk == 2 && addL5)
0779       continue;
0780     if (disk == 1 && addL6)
0781       continue;
0782     if (it < 0)
0783       disk = -disk;
0784     if (settings_.debugTracklet())
0785       edm::LogVerbatim("Tracklet") << "adding disk projection " << j << "/" << toZ_.size() << " " << disk;
0786     if (tracklet->validProj(N_LAYER + abs(disk) - 1)) {
0787       addDiskProj(tracklet, disk);
0788     }
0789   }
0790 
0791   return true;
0792 }
0793 
0794 bool TrackletCalculatorDisplaced::DDLSeeding(const Stub* innerFPGAStub,
0795                                              const L1TStub* innerStub,
0796                                              const Stub* middleFPGAStub,
0797                                              const L1TStub* middleStub,
0798                                              const Stub* outerFPGAStub,
0799                                              const L1TStub* outerStub) {
0800   if (settings_.debugTracklet())
0801     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName() << " " << layer_
0802                                  << " trying stub triplet in  (L2 D1 D2): " << innerFPGAStub->layer().value() << " "
0803                                  << middleFPGAStub->disk().value() << " " << outerFPGAStub->disk().value();
0804 
0805   int take3 = 1;  //D1D2L2
0806   unsigned ndisks = 2;
0807 
0808   double r1 = innerStub->r();
0809   double z1 = innerStub->z();
0810   double phi1 = innerStub->phi();
0811 
0812   double r2 = middleStub->r();
0813   double z2 = middleStub->z();
0814   double phi2 = middleStub->phi();
0815 
0816   double r3 = outerStub->r();
0817   double z3 = outerStub->z();
0818   double phi3 = outerStub->phi();
0819 
0820   double rinv, phi0, d0, t, z0;
0821 
0822   double phiproj[N_LAYER - 2], zproj[N_LAYER - 2], phider[N_LAYER - 2], zder[N_LAYER - 2];
0823   double phiprojdisk[N_DISK], rprojdisk[N_DISK], phiderdisk[N_DISK], rderdisk[N_DISK];
0824 
0825   exacttracklet(r1,
0826                 z1,
0827                 phi1,
0828                 r2,
0829                 z2,
0830                 phi2,
0831                 r3,
0832                 z3,
0833                 phi3,
0834                 take3,
0835                 rinv,
0836                 phi0,
0837                 d0,
0838                 t,
0839                 z0,
0840                 phiproj,
0841                 zproj,
0842                 phiprojdisk,
0843                 rprojdisk,
0844                 phider,
0845                 zder,
0846                 phiderdisk,
0847                 rderdisk);
0848 
0849   if (settings_.debugTracklet())
0850     edm::LogVerbatim("Tracklet") << __LINE__ << ":" << __FILE__ << " DLL Exact values " << innerFPGAStub->isBarrel()
0851                                  << middleFPGAStub->isBarrel() << outerFPGAStub->isBarrel() << " " << phi1 << ", " << z1
0852                                  << ", " << r1 << ", " << phi2 << ", " << z2 << ", " << r2 << ", " << phi3 << ", " << z3
0853                                  << ", " << r3 << endl;
0854 
0855   if (settings_.useapprox()) {
0856     phi1 = innerFPGAStub->phiapprox(phimin_, phimax_);
0857     z1 = innerFPGAStub->zapprox();
0858     r1 = innerFPGAStub->rapprox();
0859 
0860     phi2 = middleFPGAStub->phiapprox(phimin_, phimax_);
0861     z2 = middleFPGAStub->zapprox();
0862     r2 = middleFPGAStub->rapprox();
0863 
0864     phi3 = outerFPGAStub->phiapprox(phimin_, phimax_);
0865     z3 = outerFPGAStub->zapprox();
0866     r3 = outerFPGAStub->rapprox();
0867   }
0868 
0869   if (settings_.debugTracklet())
0870     edm::LogVerbatim("Tracklet") << __LINE__ << ":" << __FILE__ << "DLL Approx values " << innerFPGAStub->isBarrel()
0871                                  << middleFPGAStub->isBarrel() << outerFPGAStub->isBarrel() << " " << phi1 << ", " << z1
0872                                  << ", " << r1 << ", " << phi2 << ", " << z2 << ", " << r2 << ", " << phi3 << ", " << z3
0873                                  << ", " << r3 << endl;
0874 
0875   double rinvapprox, phi0approx, d0approx, tapprox, z0approx;
0876   double phiprojapprox[N_LAYER - 2], zprojapprox[N_LAYER - 2], phiderapprox[N_LAYER - 2], zderapprox[N_LAYER - 2];
0877   double phiprojdiskapprox[N_DISK], rprojdiskapprox[N_DISK];
0878   double phiderdiskapprox[N_DISK], rderdiskapprox[N_DISK];
0879 
0880   //TODO: implement the actual integer calculation
0881   if (settings_.useapprox()) {
0882     approxtracklet(r1,
0883                    z1,
0884                    phi1,
0885                    r2,
0886                    z2,
0887                    phi2,
0888                    r3,
0889                    z3,
0890                    phi3,
0891                    take3,
0892                    ndisks,
0893                    rinvapprox,
0894                    phi0approx,
0895                    d0approx,
0896                    tapprox,
0897                    z0approx,
0898                    phiprojapprox,
0899                    zprojapprox,
0900                    phiderapprox,
0901                    zderapprox,
0902                    phiprojdiskapprox,
0903                    rprojdiskapprox,
0904                    phiderdiskapprox,
0905                    rderdiskapprox);
0906   } else {
0907     rinvapprox = rinv;
0908     phi0approx = phi0;
0909     d0approx = d0;
0910     tapprox = t;
0911     z0approx = z0;
0912 
0913     for (unsigned int i = 0; i < toR_.size(); ++i) {
0914       phiprojapprox[i] = phiproj[i];
0915       zprojapprox[i] = zproj[i];
0916       phiderapprox[i] = phider[i];
0917       zderapprox[i] = zder[i];
0918     }
0919 
0920     for (unsigned int i = 0; i < toZ_.size(); ++i) {
0921       phiprojdiskapprox[i] = phiprojdisk[i];
0922       rprojdiskapprox[i] = rprojdisk[i];
0923       phiderdiskapprox[i] = phiderdisk[i];
0924       rderdiskapprox[i] = rderdisk[i];
0925     }
0926   }
0927 
0928   //store the approcximate results
0929   if (settings_.debugTracklet()) {
0930     edm::LogVerbatim("Tracklet") << "rinvapprox: " << rinvapprox << " rinv: " << rinv << endl;
0931     edm::LogVerbatim("Tracklet") << "phi0approx: " << phi0approx << " phi0: " << phi0 << endl;
0932     edm::LogVerbatim("Tracklet") << "d0approx: " << d0approx << " d0: " << d0 << endl;
0933     edm::LogVerbatim("Tracklet") << "tapprox: " << tapprox << " t: " << t << endl;
0934     edm::LogVerbatim("Tracklet") << "z0approx: " << z0approx << " z0: " << z0 << endl;
0935   }
0936 
0937   for (unsigned int i = 0; i < toR_.size(); ++i) {
0938     if (settings_.debugTracklet()) {
0939       edm::LogVerbatim("Tracklet") << "phiprojapprox[" << i << "]: " << phiprojapprox[i] << " phiproj[" << i
0940                                    << "]: " << phiproj[i] << endl;
0941       edm::LogVerbatim("Tracklet") << "zprojapprox[" << i << "]: " << zprojapprox[i] << " zproj[" << i
0942                                    << "]: " << zproj[i] << endl;
0943       edm::LogVerbatim("Tracklet") << "phiderapprox[" << i << "]: " << phiderapprox[i] << " phider[" << i
0944                                    << "]: " << phider[i] << endl;
0945       edm::LogVerbatim("Tracklet") << "zderapprox[" << i << "]: " << zderapprox[i] << " zder[" << i << "]: " << zder[i]
0946                                    << endl;
0947     }
0948   }
0949 
0950   for (unsigned int i = 0; i < toZ_.size(); ++i) {
0951     if (settings_.debugTracklet()) {
0952       edm::LogVerbatim("Tracklet") << "phiprojdiskapprox[" << i << "]: " << phiprojdiskapprox[i] << " phiprojdisk[" << i
0953                                    << "]: " << phiprojdisk[i] << endl;
0954       edm::LogVerbatim("Tracklet") << "rprojdiskapprox[" << i << "]: " << rprojdiskapprox[i] << " rprojdisk[" << i
0955                                    << "]: " << rprojdisk[i] << endl;
0956       edm::LogVerbatim("Tracklet") << "phiderdiskapprox[" << i << "]: " << phiderdiskapprox[i] << " phiderdisk[" << i
0957                                    << "]: " << phiderdisk[i] << endl;
0958       edm::LogVerbatim("Tracklet") << "rderdiskapprox[" << i << "]: " << rderdiskapprox[i] << " rderdisk[" << i
0959                                    << "]: " << rderdisk[i] << endl;
0960     }
0961   }
0962 
0963   //now binary
0964   double krinv = settings_.kphi1() / settings_.kr() * pow(2, settings_.rinv_shift()),
0965          kphi0 = settings_.kphi1() * pow(2, settings_.phi0_shift()),
0966          kt = settings_.kz() / settings_.kr() * pow(2, settings_.t_shift()),
0967          kz0 = settings_.kz() * pow(2, settings_.z0_shift()),
0968          kphiproj = settings_.kphi1() * pow(2, settings_.SS_phiL_shift()),
0969          kphider = settings_.kphi1() / settings_.kr() * pow(2, settings_.SS_phiderL_shift()),
0970          kzproj = settings_.kz() * pow(2, settings_.PS_zL_shift()),
0971          kzder = settings_.kz() / settings_.kr() * pow(2, settings_.PS_zderL_shift()),
0972          kphiprojdisk = settings_.kphi1() * pow(2, settings_.SS_phiD_shift()),
0973          kphiderdisk = settings_.kphi1() / settings_.kr() * pow(2, settings_.SS_phiderD_shift()),
0974          krprojdisk = settings_.kr() * pow(2, settings_.PS_rD_shift()),
0975          krderdisk = settings_.kr() / settings_.kz() * pow(2, settings_.PS_rderD_shift());
0976 
0977   int irinv, iphi0, id0, it, iz0;
0978   int iphiproj[N_LAYER - 2], izproj[N_LAYER - 2], iphider[N_LAYER - 2], izder[N_LAYER - 2];
0979   int iphiprojdisk[N_DISK], irprojdisk[N_DISK], iphiderdisk[N_DISK], irderdisk[N_DISK];
0980 
0981   //store the binary results
0982   irinv = rinvapprox / krinv;
0983   iphi0 = phi0approx / kphi0;
0984   id0 = d0approx / settings_.kd0();
0985   it = tapprox / kt;
0986   iz0 = z0approx / kz0;
0987 
0988   bool success = true;
0989   if (std::abs(rinvapprox) > settings_.rinvcut()) {
0990     if (settings_.debugTracklet())
0991       edm::LogVerbatim("Tracklet") << "TrackletCalculator::DDL Seeding irinv too large: " << rinvapprox << "(" << irinv
0992                                    << ")";
0993     success = false;
0994   }
0995   if (std::abs(z0approx) > settings_.disp_z0cut()) {
0996     if (settings_.debugTracklet())
0997       edm::LogVerbatim("Tracklet") << "Failed tracklet z0 cut " << z0approx;
0998     success = false;
0999   }
1000   if (std::abs(d0approx) > settings_.maxd0()) {
1001     if (settings_.debugTracklet())
1002       edm::LogVerbatim("Tracklet") << "Failed tracklet d0 cut " << d0approx;
1003     success = false;
1004   }
1005 
1006   if (!success)
1007     return false;
1008 
1009   double phicritapprox = phi0approx - asin(0.5 * settings_.rcrit() * rinvapprox);
1010   int phicrit = iphi0 - 2 * irinv;
1011 
1012   int iphicritmincut = settings_.phicritminmc() / globals_->ITC_L1L2()->phi0_final.K();
1013   int iphicritmaxcut = settings_.phicritmaxmc() / globals_->ITC_L1L2()->phi0_final.K();
1014 
1015   bool keepapprox = (phicritapprox > settings_.phicritminmc()) && (phicritapprox < settings_.phicritmaxmc()),
1016        keep = (phicrit > iphicritmincut) && (phicrit < iphicritmaxcut);
1017 
1018   if (settings_.debugTracklet())
1019     if (keep && !keepapprox)
1020       edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced::DDLSeeding tracklet kept with exact phicrit cut "
1021                                       "but not approximate, phicritapprox: "
1022                                    << phicritapprox;
1023   if (settings_.usephicritapprox()) {
1024     if (!keepapprox)
1025       return false;
1026   } else {
1027     if (!keep)
1028       return false;
1029   }
1030 
1031   Projection projs[N_LAYER + N_DISK];
1032 
1033   for (unsigned int i = 0; i < toR_.size(); ++i) {
1034     iphiproj[i] = phiprojapprox[i] / kphiproj;
1035     izproj[i] = zprojapprox[i] / kzproj;
1036 
1037     iphider[i] = phiderapprox[i] / kphider;
1038     izder[i] = zderapprox[i] / kzder;
1039 
1040     //check that z projection in range
1041     if (izproj[i] < -(1 << (settings_.nzbitsstub(0) - 1)))
1042       continue;
1043     if (izproj[i] >= (1 << (settings_.nzbitsstub(0) - 1)))
1044       continue;
1045 
1046     //check that phi projection in range
1047     if (iphiproj[i] >= (1 << settings_.nphibitsstub(5)) - 1)
1048       continue;
1049     if (iphiproj[i] <= 0)
1050       continue;
1051 
1052     if (rproj_[i] < settings_.rPS2S()) {
1053       iphiproj[i] >>= (settings_.nphibitsstub(5) - settings_.nphibitsstub(0));
1054     } else {
1055       izproj[i] >>= (settings_.nzbitsstub(0) - settings_.nzbitsstub(5));
1056     }
1057 
1058     if (rproj_[i] < settings_.rPS2S()) {
1059       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL123() - 1)))
1060         iphider[i] = -(1 << (settings_.nbitsphiprojderL123() - 1));
1061       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL123() - 1)))
1062         iphider[i] = (1 << (settings_.nbitsphiprojderL123() - 1)) - 1;
1063     } else {
1064       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL456() - 1)))
1065         iphider[i] = -(1 << (settings_.nbitsphiprojderL456() - 1));
1066       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL456() - 1)))
1067         iphider[i] = (1 << (settings_.nbitsphiprojderL456() - 1)) - 1;
1068     }
1069 
1070     projs[lproj_[i] - 1].init(settings_,
1071                               lproj_[i] - 1,
1072                               iphiproj[i],
1073                               izproj[i],
1074                               iphider[i],
1075                               izder[i],
1076                               phiproj[i],
1077                               zproj[i],
1078                               phider[i],
1079                               zder[i],
1080                               phiprojapprox[i],
1081                               zprojapprox[i],
1082                               phiderapprox[i],
1083                               zderapprox[i],
1084                               false);
1085   }
1086 
1087   if (std::abs(it * kt) > 1.0) {
1088     for (unsigned int i = 0; i < toZ_.size(); ++i) {
1089       iphiprojdisk[i] = phiprojdiskapprox[i] / kphiprojdisk;
1090       irprojdisk[i] = rprojdiskapprox[i] / krprojdisk;
1091 
1092       iphiderdisk[i] = phiderdiskapprox[i] / kphiderdisk;
1093       irderdisk[i] = rderdiskapprox[i] / krderdisk;
1094 
1095       if (iphiprojdisk[i] <= 0)
1096         continue;
1097       if (iphiprojdisk[i] >= (1 << settings_.nphibitsstub(0)) - 1)
1098         continue;
1099 
1100       if (irprojdisk[i] < settings_.rmindisk() / krprojdisk || irprojdisk[i] > settings_.rmaxdisk() / krprojdisk)
1101         continue;
1102 
1103       projs[N_LAYER + i + 2].init(settings_,
1104                                   N_LAYER + i + 2,
1105                                   iphiprojdisk[i],
1106                                   irprojdisk[i],
1107                                   iphiderdisk[i],
1108                                   irderdisk[i],
1109                                   phiprojdisk[i],
1110                                   rprojdisk[i],
1111                                   phiderdisk[i],
1112                                   rderdisk[i],
1113                                   phiprojdiskapprox[i],
1114                                   rprojdiskapprox[i],
1115                                   phiderdisk[i],
1116                                   rderdisk[i],
1117                                   false);
1118     }
1119   }
1120 
1121   if (settings_.writeMonitorData("TrackletPars")) {
1122     globals_->ofstream("trackletpars.txt")
1123         << layer_ << " , " << rinv << " , " << rinvapprox << " , " << phi0 << " , " << phi0approx << " , " << t << " , "
1124         << tapprox << " , " << z0 << " , " << z0approx << " , " << d0 << " , " << d0approx << endl;
1125   }
1126 
1127   Tracklet* tracklet = new Tracklet(settings_,
1128                                     iSeed_,
1129                                     innerFPGAStub,
1130                                     middleFPGAStub,
1131                                     outerFPGAStub,
1132                                     rinv,
1133                                     phi0,
1134                                     d0,
1135                                     z0,
1136                                     t,
1137                                     rinvapprox,
1138                                     phi0approx,
1139                                     d0approx,
1140                                     z0approx,
1141                                     tapprox,
1142                                     irinv,
1143                                     iphi0,
1144                                     id0,
1145                                     iz0,
1146                                     it,
1147                                     projs,
1148                                     true);
1149 
1150   if (settings_.debugTracklet())
1151     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName()
1152                                  << " Found DDL tracklet in sector = " << iSector_ << " phi0 = " << phi0;
1153 
1154   tracklet->setTrackletIndex(trackletpars_->nTracklets());
1155   tracklet->setTCIndex(TCIndex_);
1156 
1157   if (settings_.writeMonitorData("Seeds")) {
1158     ofstream fout("seeds.txt", ofstream::app);
1159     fout << __FILE__ << ":" << __LINE__ << " " << name_ << "_" << iSector_ << " " << tracklet->getISeed() << endl;
1160     fout.close();
1161   }
1162   trackletpars_->addTracklet(tracklet);
1163 
1164   for (unsigned int j = 0; j < toR_.size(); j++) {
1165     if (settings_.debugTracklet())
1166       edm::LogVerbatim("Tracklet") << "adding layer projection " << j << "/" << toR_.size() << " " << lproj_[j] << " "
1167                                    << tracklet->validProj(lproj_[j] - 1);
1168     if (tracklet->validProj(lproj_[j] - 1)) {
1169       addLayerProj(tracklet, lproj_[j]);
1170     }
1171   }
1172 
1173   for (unsigned int j = 0; j < toZ_.size(); j++) {
1174     int disk = dproj_[j];
1175     if (disk == 0)
1176       continue;
1177     if (it < 0)
1178       disk = -disk;
1179     if (settings_.debugTracklet())
1180       edm::LogVerbatim("Tracklet") << "adding disk projection " << j << "/" << toZ_.size() << " " << disk << " "
1181                                    << tracklet->validProj(N_LAYER + abs(disk) - 1);
1182     if (tracklet->validProj(N_LAYER + abs(disk) - 1)) {
1183       addDiskProj(tracklet, disk);
1184     }
1185   }
1186 
1187   return true;
1188 }
1189 
1190 bool TrackletCalculatorDisplaced::LLDSeeding(const Stub* innerFPGAStub,
1191                                              const L1TStub* innerStub,
1192                                              const Stub* middleFPGAStub,
1193                                              const L1TStub* middleStub,
1194                                              const Stub* outerFPGAStub,
1195                                              const L1TStub* outerStub) {
1196   if (settings_.debugTracklet())
1197     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName() << " " << layer_
1198                                  << " trying stub triplet in  (L2L3D1): " << middleFPGAStub->layer().value() << " "
1199                                  << outerFPGAStub->layer().value() << " " << innerFPGAStub->disk().value();
1200 
1201   int take3 = 0;  //L2L3D1
1202   unsigned ndisks = 1;
1203 
1204   double r3 = innerStub->r();
1205   double z3 = innerStub->z();
1206   double phi3 = innerStub->phi();
1207 
1208   double r1 = middleStub->r();
1209   double z1 = middleStub->z();
1210   double phi1 = middleStub->phi();
1211 
1212   double r2 = outerStub->r();
1213   double z2 = outerStub->z();
1214   double phi2 = outerStub->phi();
1215 
1216   double rinv, phi0, d0, t, z0;
1217 
1218   double phiproj[N_LAYER - 2], zproj[N_LAYER - 2], phider[N_LAYER - 2], zder[N_LAYER - 2];
1219   double phiprojdisk[N_DISK], rprojdisk[N_DISK], phiderdisk[N_DISK], rderdisk[N_DISK];
1220 
1221   exacttracklet(r1,
1222                 z1,
1223                 phi1,
1224                 r2,
1225                 z2,
1226                 phi2,
1227                 r3,
1228                 z3,
1229                 phi3,
1230                 take3,
1231                 rinv,
1232                 phi0,
1233                 d0,
1234                 t,
1235                 z0,
1236                 phiproj,
1237                 zproj,
1238                 phiprojdisk,
1239                 rprojdisk,
1240                 phider,
1241                 zder,
1242                 phiderdisk,
1243                 rderdisk);
1244 
1245   if (settings_.debugTracklet())
1246     edm::LogVerbatim("Tracklet") << __LINE__ << ":" << __FILE__ << "LLD Exact values " << innerFPGAStub->isBarrel()
1247                                  << middleFPGAStub->isBarrel() << outerFPGAStub->isBarrel() << " " << phi3 << ", " << z3
1248                                  << ", " << r3 << ", " << phi1 << ", " << z1 << ", " << r1 << ", " << phi2 << ", " << z2
1249                                  << ", " << r2 << endl;
1250 
1251   if (settings_.useapprox()) {
1252     phi3 = innerFPGAStub->phiapprox(phimin_, phimax_);
1253     z3 = innerFPGAStub->zapprox();
1254     r3 = innerFPGAStub->rapprox();
1255 
1256     phi1 = middleFPGAStub->phiapprox(phimin_, phimax_);
1257     z1 = middleFPGAStub->zapprox();
1258     r1 = middleFPGAStub->rapprox();
1259 
1260     phi2 = outerFPGAStub->phiapprox(phimin_, phimax_);
1261     z2 = outerFPGAStub->zapprox();
1262     r2 = outerFPGAStub->rapprox();
1263   }
1264 
1265   if (settings_.debugTracklet())
1266     edm::LogVerbatim("Tracklet") << __LINE__ << ":" << __FILE__ << "LLD approx values " << innerFPGAStub->isBarrel()
1267                                  << middleFPGAStub->isBarrel() << outerFPGAStub->isBarrel() << " " << phi3 << ", " << z3
1268                                  << ", " << r3 << ", " << phi1 << ", " << z1 << ", " << r1 << ", " << phi2 << ", " << z2
1269                                  << ", " << r2 << endl;
1270 
1271   double rinvapprox, phi0approx, d0approx, tapprox, z0approx;
1272   double phiprojapprox[N_LAYER - 2], zprojapprox[N_LAYER - 2], phiderapprox[N_LAYER - 2], zderapprox[N_LAYER - 2];
1273   double phiprojdiskapprox[N_DISK], rprojdiskapprox[N_DISK];
1274   double phiderdiskapprox[N_DISK], rderdiskapprox[N_DISK];
1275 
1276   //TODO: implement the actual integer calculation
1277   if (settings_.useapprox()) {
1278     approxtracklet(r1,
1279                    z1,
1280                    phi1,
1281                    r2,
1282                    z2,
1283                    phi2,
1284                    r3,
1285                    z3,
1286                    phi3,
1287                    take3,
1288                    ndisks,
1289                    rinvapprox,
1290                    phi0approx,
1291                    d0approx,
1292                    tapprox,
1293                    z0approx,
1294                    phiprojapprox,
1295                    zprojapprox,
1296                    phiderapprox,
1297                    zderapprox,
1298                    phiprojdiskapprox,
1299                    rprojdiskapprox,
1300                    phiderdiskapprox,
1301                    rderdiskapprox);
1302   } else {
1303     rinvapprox = rinv;
1304     phi0approx = phi0;
1305     d0approx = d0;
1306     tapprox = t;
1307     z0approx = z0;
1308 
1309     for (unsigned int i = 0; i < toR_.size(); ++i) {
1310       phiprojapprox[i] = phiproj[i];
1311       zprojapprox[i] = zproj[i];
1312       phiderapprox[i] = phider[i];
1313       zderapprox[i] = zder[i];
1314     }
1315 
1316     for (unsigned int i = 0; i < toZ_.size(); ++i) {
1317       phiprojdiskapprox[i] = phiprojdisk[i];
1318       rprojdiskapprox[i] = rprojdisk[i];
1319       phiderdiskapprox[i] = phiderdisk[i];
1320       rderdiskapprox[i] = rderdisk[i];
1321     }
1322   }
1323 
1324   //store the approcximate results
1325   if (settings_.debugTracklet()) {
1326     edm::LogVerbatim("Tracklet") << "rinvapprox: " << rinvapprox << " rinv: " << rinv << endl;
1327     edm::LogVerbatim("Tracklet") << "phi0approx: " << phi0approx << " phi0: " << phi0 << endl;
1328     edm::LogVerbatim("Tracklet") << "d0approx: " << d0approx << " d0: " << d0 << endl;
1329     edm::LogVerbatim("Tracklet") << "tapprox: " << tapprox << " t: " << t << endl;
1330     edm::LogVerbatim("Tracklet") << "z0approx: " << z0approx << " z0: " << z0 << endl;
1331   }
1332 
1333   for (unsigned int i = 0; i < toR_.size(); ++i) {
1334     if (settings_.debugTracklet()) {
1335       edm::LogVerbatim("Tracklet") << "phiprojapprox[" << i << "]: " << phiprojapprox[i] << " phiproj[" << i
1336                                    << "]: " << phiproj[i] << endl;
1337       edm::LogVerbatim("Tracklet") << "zprojapprox[" << i << "]: " << zprojapprox[i] << " zproj[" << i
1338                                    << "]: " << zproj[i] << endl;
1339       edm::LogVerbatim("Tracklet") << "phiderapprox[" << i << "]: " << phiderapprox[i] << " phider[" << i
1340                                    << "]: " << phider[i] << endl;
1341       edm::LogVerbatim("Tracklet") << "zderapprox[" << i << "]: " << zderapprox[i] << " zder[" << i << "]: " << zder[i]
1342                                    << endl;
1343     }
1344   }
1345 
1346   for (unsigned int i = 0; i < toZ_.size(); ++i) {
1347     if (settings_.debugTracklet()) {
1348       edm::LogVerbatim("Tracklet") << "phiprojdiskapprox[" << i << "]: " << phiprojdiskapprox[i] << " phiprojdisk[" << i
1349                                    << "]: " << phiprojdisk[i] << endl;
1350       edm::LogVerbatim("Tracklet") << "rprojdiskapprox[" << i << "]: " << rprojdiskapprox[i] << " rprojdisk[" << i
1351                                    << "]: " << rprojdisk[i] << endl;
1352       edm::LogVerbatim("Tracklet") << "phiderdiskapprox[" << i << "]: " << phiderdiskapprox[i] << " phiderdisk[" << i
1353                                    << "]: " << phiderdisk[i] << endl;
1354       edm::LogVerbatim("Tracklet") << "rderdiskapprox[" << i << "]: " << rderdiskapprox[i] << " rderdisk[" << i
1355                                    << "]: " << rderdisk[i] << endl;
1356     }
1357   }
1358 
1359   //now binary
1360   double krinv = settings_.kphi1() / settings_.kr() * pow(2, settings_.rinv_shift()),
1361          kphi0 = settings_.kphi1() * pow(2, settings_.phi0_shift()),
1362          kt = settings_.kz() / settings_.kr() * pow(2, settings_.t_shift()),
1363          kz0 = settings_.kz() * pow(2, settings_.z0_shift()),
1364          kphiproj = settings_.kphi1() * pow(2, settings_.SS_phiL_shift()),
1365          kphider = settings_.kphi1() / settings_.kr() * pow(2, settings_.SS_phiderL_shift()),
1366          kzproj = settings_.kz() * pow(2, settings_.PS_zL_shift()),
1367          kzder = settings_.kz() / settings_.kr() * pow(2, settings_.PS_zderL_shift()),
1368          kphiprojdisk = settings_.kphi1() * pow(2, settings_.SS_phiD_shift()),
1369          kphiderdisk = settings_.kphi1() / settings_.kr() * pow(2, settings_.SS_phiderD_shift()),
1370          krprojdisk = settings_.kr() * pow(2, settings_.PS_rD_shift()),
1371          krderdisk = settings_.kr() / settings_.kz() * pow(2, settings_.PS_rderD_shift());
1372 
1373   int irinv, iphi0, id0, it, iz0;
1374   int iphiproj[N_LAYER - 2], izproj[N_LAYER - 2], iphider[N_LAYER - 2], izder[N_LAYER - 2];
1375   int iphiprojdisk[N_DISK], irprojdisk[N_DISK], iphiderdisk[N_DISK], irderdisk[N_DISK];
1376 
1377   //store the binary results
1378   irinv = rinvapprox / krinv;
1379   iphi0 = phi0approx / kphi0;
1380   id0 = d0approx / settings_.kd0();
1381   it = tapprox / kt;
1382   iz0 = z0approx / kz0;
1383 
1384   bool success = true;
1385   if (std::abs(rinvapprox) > settings_.rinvcut()) {
1386     if (settings_.debugTracklet())
1387       edm::LogVerbatim("Tracklet") << "TrackletCalculator:: LLD Seeding irinv too large: " << rinvapprox << "(" << irinv
1388                                    << ")";
1389     success = false;
1390   }
1391   if (std::abs(z0approx) > settings_.disp_z0cut()) {
1392     if (settings_.debugTracklet())
1393       edm::LogVerbatim("Tracklet") << "Failed tracklet z0 cut " << z0approx;
1394     success = false;
1395   }
1396   if (std::abs(d0approx) > settings_.maxd0()) {
1397     if (settings_.debugTracklet())
1398       edm::LogVerbatim("Tracklet") << "Failed tracklet d0 cut " << d0approx;
1399     success = false;
1400   }
1401 
1402   if (!success)
1403     return false;
1404 
1405   double phicritapprox = phi0approx - asin(0.5 * settings_.rcrit() * rinvapprox);
1406   int phicrit = iphi0 - 2 * irinv;
1407 
1408   int iphicritmincut = settings_.phicritminmc() / globals_->ITC_L1L2()->phi0_final.K();
1409   int iphicritmaxcut = settings_.phicritmaxmc() / globals_->ITC_L1L2()->phi0_final.K();
1410 
1411   bool keepapprox = (phicritapprox > settings_.phicritminmc()) && (phicritapprox < settings_.phicritmaxmc()),
1412        keep = (phicrit > iphicritmincut) && (phicrit < iphicritmaxcut);
1413 
1414   if (settings_.debugTracklet())
1415     if (keep && !keepapprox)
1416       edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced::LLDSeeding tracklet kept with exact phicrit cut "
1417                                       "but not approximate, phicritapprox: "
1418                                    << phicritapprox;
1419   if (settings_.usephicritapprox()) {
1420     if (!keepapprox)
1421       return false;
1422   } else {
1423     if (!keep)
1424       return false;
1425   }
1426 
1427   Projection projs[N_LAYER + N_DISK];
1428 
1429   for (unsigned int i = 0; i < toR_.size(); ++i) {
1430     iphiproj[i] = phiprojapprox[i] / kphiproj;
1431     izproj[i] = zprojapprox[i] / kzproj;
1432 
1433     iphider[i] = phiderapprox[i] / kphider;
1434     izder[i] = zderapprox[i] / kzder;
1435 
1436     if (izproj[i] < -(1 << (settings_.nzbitsstub(0) - 1)))
1437       continue;
1438     if (izproj[i] >= (1 << (settings_.nzbitsstub(0) - 1)))
1439       continue;
1440 
1441     //this is left from the original....
1442     if (iphiproj[i] >= (1 << settings_.nphibitsstub(5)) - 1)
1443       continue;
1444     if (iphiproj[i] <= 0)
1445       continue;
1446 
1447     if (rproj_[i] < settings_.rPS2S()) {
1448       iphiproj[i] >>= (settings_.nphibitsstub(5) - settings_.nphibitsstub(0));
1449     } else {
1450       izproj[i] >>= (settings_.nzbitsstub(0) - settings_.nzbitsstub(5));
1451     }
1452 
1453     if (rproj_[i] < settings_.rPS2S()) {
1454       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL123() - 1)))
1455         iphider[i] = -(1 << (settings_.nbitsphiprojderL123() - 1));
1456       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL123() - 1)))
1457         iphider[i] = (1 << (settings_.nbitsphiprojderL123() - 1)) - 1;
1458     } else {
1459       if (iphider[i] < -(1 << (settings_.nbitsphiprojderL456() - 1)))
1460         iphider[i] = -(1 << (settings_.nbitsphiprojderL456() - 1));
1461       if (iphider[i] >= (1 << (settings_.nbitsphiprojderL456() - 1)))
1462         iphider[i] = (1 << (settings_.nbitsphiprojderL456() - 1)) - 1;
1463     }
1464 
1465     projs[lproj_[i] - 1].init(settings_,
1466                               lproj_[i] - 1,
1467                               iphiproj[i],
1468                               izproj[i],
1469                               iphider[i],
1470                               izder[i],
1471                               phiproj[i],
1472                               zproj[i],
1473                               phider[i],
1474                               zder[i],
1475                               phiprojapprox[i],
1476                               zprojapprox[i],
1477                               phiderapprox[i],
1478                               zderapprox[i],
1479                               false);
1480   }
1481 
1482   if (std::abs(it * kt) > 1.0) {
1483     for (unsigned int i = 0; i < toZ_.size(); ++i) {
1484       iphiprojdisk[i] = phiprojdiskapprox[i] / kphiprojdisk;
1485       irprojdisk[i] = rprojdiskapprox[i] / krprojdisk;
1486 
1487       iphiderdisk[i] = phiderdiskapprox[i] / kphiderdisk;
1488       irderdisk[i] = rderdiskapprox[i] / krderdisk;
1489 
1490       //Check phi range of projection
1491       if (iphiprojdisk[i] <= 0)
1492         continue;
1493       if (iphiprojdisk[i] >= (1 << settings_.nphibitsstub(0)) - 1)
1494         continue;
1495 
1496       //Check r range of projection
1497       if (irprojdisk[i] < settings_.rmindisk() / krprojdisk || irprojdisk[i] > settings_.rmaxdisk() / krprojdisk)
1498         continue;
1499 
1500       projs[N_LAYER + i + 1].init(settings_,
1501                                   N_LAYER + i + 1,
1502                                   iphiprojdisk[i],
1503                                   irprojdisk[i],
1504                                   iphiderdisk[i],
1505                                   irderdisk[i],
1506                                   phiprojdisk[i],
1507                                   rprojdisk[i],
1508                                   phiderdisk[i],
1509                                   rderdisk[i],
1510                                   phiprojdiskapprox[i],
1511                                   rprojdiskapprox[i],
1512                                   phiderdisk[i],
1513                                   rderdisk[i],
1514                                   false);
1515     }
1516   }
1517 
1518   if (settings_.writeMonitorData("TrackletPars")) {
1519     globals_->ofstream("trackletpars.txt")
1520         << layer_ << " , " << rinv << " , " << rinvapprox << " , " << phi0 << " , " << phi0approx << " , " << t << " , "
1521         << tapprox << " , " << z0 << " , " << z0approx << " , " << d0 << " , " << d0approx << endl;
1522   }
1523 
1524   Tracklet* tracklet = new Tracklet(settings_,
1525                                     iSeed_,
1526                                     innerFPGAStub,
1527                                     middleFPGAStub,
1528                                     outerFPGAStub,
1529                                     rinv,
1530                                     phi0,
1531                                     d0,
1532                                     z0,
1533                                     t,
1534                                     rinvapprox,
1535                                     phi0approx,
1536                                     d0approx,
1537                                     z0approx,
1538                                     tapprox,
1539                                     irinv,
1540                                     iphi0,
1541                                     id0,
1542                                     iz0,
1543                                     it,
1544                                     projs,
1545                                     false);
1546 
1547   if (settings_.debugTracklet())
1548     edm::LogVerbatim("Tracklet") << "TrackletCalculatorDisplaced " << getName()
1549                                  << " Found LLD tracklet in sector = " << iSector_ << " phi0 = " << phi0;
1550 
1551   tracklet->setTrackletIndex(trackletpars_->nTracklets());
1552   tracklet->setTCIndex(TCIndex_);
1553 
1554   if (settings_.writeMonitorData("Seeds")) {
1555     ofstream fout("seeds.txt", ofstream::app);
1556     fout << __FILE__ << ":" << __LINE__ << " " << name_ << "_" << iSector_ << " " << tracklet->getISeed() << endl;
1557     fout.close();
1558   }
1559   trackletpars_->addTracklet(tracklet);
1560 
1561   for (unsigned int j = 0; j < toR_.size(); j++) {
1562     if (settings_.debugTracklet())
1563       edm::LogVerbatim("Tracklet") << "adding layer projection " << j << "/" << toR_.size() << " " << lproj_[j];
1564     if (tracklet->validProj(lproj_[j] - 1)) {
1565       addLayerProj(tracklet, lproj_[j]);
1566     }
1567   }
1568 
1569   for (unsigned int j = 0; j < toZ_.size(); j++) {
1570     int disk = dproj_[j];
1571     if (disk == 0)
1572       continue;
1573     if (it < 0)
1574       disk = -disk;
1575     if (settings_.debugTracklet())
1576       edm::LogVerbatim("Tracklet") << "adding disk projection " << j << "/" << toZ_.size() << " " << disk;
1577     if (tracklet->validProj(N_LAYER + abs(disk) - 1)) {
1578       addDiskProj(tracklet, disk);
1579     }
1580   }
1581 
1582   return true;
1583 }
1584 
1585 void TrackletCalculatorDisplaced::exactproj(double rproj,
1586                                             double rinv,
1587                                             double phi0,
1588                                             double d0,
1589                                             double t,
1590                                             double z0,
1591                                             double r0,
1592                                             double& phiproj,
1593                                             double& zproj,
1594                                             double& phider,
1595                                             double& zder) {
1596   double rho = 1 / rinv;
1597   if (rho < 0) {
1598     r0 = -r0;
1599   }
1600   phiproj = phi0 - asin((rproj * rproj + r0 * r0 - rho * rho) / (2 * rproj * r0));
1601   double beta = acos((rho * rho + r0 * r0 - rproj * rproj) / (2 * r0 * rho));
1602   zproj = z0 + t * std::abs(rho * beta);
1603 
1604   //not exact, but close
1605   phider = -0.5 * rinv / sqrt(1 - pow(0.5 * rproj * rinv, 2)) + d0 / (rproj * rproj);
1606   zder = t / sqrt(1 - pow(0.5 * rproj * rinv, 2));
1607 
1608   if (settings_.debugTracklet())
1609     edm::LogVerbatim("Tracklet") << "exact proj layer at " << rproj << " : " << phiproj << " " << zproj;
1610 }
1611 
1612 void TrackletCalculatorDisplaced::exactprojdisk(double zproj,
1613                                                 double rinv,
1614                                                 double,
1615                                                 double,  //phi0 and d0 are not used.
1616                                                 double t,
1617                                                 double z0,
1618                                                 double x0,
1619                                                 double y0,
1620                                                 double& phiproj,
1621                                                 double& rproj,
1622                                                 double& phider,
1623                                                 double& rder) {
1624   //protect against t=0
1625   if (std::abs(t) < 0.1)
1626     t = 0.1;
1627   if (t < 0)
1628     zproj = -zproj;
1629   double rho = std::abs(1 / rinv);
1630   double beta = (zproj - z0) / (t * rho);
1631   double phiV = atan2(-y0, -x0);
1632   double c = rinv > 0 ? -1 : 1;
1633 
1634   double x = x0 + rho * cos(phiV + c * beta);
1635   double y = y0 + rho * sin(phiV + c * beta);
1636 
1637   phiproj = atan2(y, x);
1638 
1639   phiproj = reco::reduceRange(phiproj - phimin_);
1640 
1641   rproj = sqrt(x * x + y * y);
1642 
1643   phider = c / t / (x * x + y * y) * (rho + x0 * cos(phiV + c * beta) + y0 * sin(phiV + c * beta));
1644   rder = c / t / rproj * (y0 * cos(phiV + c * beta) - x0 * sin(phiV + c * beta));
1645 
1646   if (settings_.debugTracklet())
1647     edm::LogVerbatim("Tracklet") << "exact proj disk at" << zproj << " : " << phiproj << " " << rproj;
1648 }
1649 
1650 void TrackletCalculatorDisplaced::exacttracklet(double r1,
1651                                                 double z1,
1652                                                 double phi1,
1653                                                 double r2,
1654                                                 double z2,
1655                                                 double phi2,
1656                                                 double r3,
1657                                                 double z3,
1658                                                 double phi3,
1659                                                 int take3,
1660                                                 double& rinv,
1661                                                 double& phi0,
1662                                                 double& d0,
1663                                                 double& t,
1664                                                 double& z0,
1665                                                 double phiproj[N_DISK],
1666                                                 double zproj[N_DISK],
1667                                                 double phiprojdisk[N_DISK],
1668                                                 double rprojdisk[N_DISK],
1669                                                 double phider[N_DISK],
1670                                                 double zder[N_DISK],
1671                                                 double phiderdisk[N_DISK],
1672                                                 double rderdisk[N_DISK]) {
1673   //two lines perpendicular to the 1->2 and 2->3
1674   double x1 = r1 * cos(phi1);
1675   double x2 = r2 * cos(phi2);
1676   double x3 = r3 * cos(phi3);
1677 
1678   double y1 = r1 * sin(phi1);
1679   double y2 = r2 * sin(phi2);
1680   double y3 = r3 * sin(phi3);
1681 
1682   double dy21 = y2 - y1;
1683   double dy32 = y3 - y2;
1684 
1685   //Hack to protect against dividing by zero
1686   //code should be rewritten to avoid this
1687   if (dy21 == 0.0)
1688     dy21 = 1e-9;
1689   if (dy32 == 0.0)
1690     dy32 = 1e-9;
1691 
1692   double k1 = -(x2 - x1) / dy21;
1693   double k2 = -(x3 - x2) / dy32;
1694   double b1 = 0.5 * (y2 + y1) - 0.5 * (x1 + x2) * k1;
1695   double b2 = 0.5 * (y3 + y2) - 0.5 * (x2 + x3) * k2;
1696   //their intersection gives the center of the circle
1697   double y0 = (b1 * k2 - b2 * k1) / (k2 - k1);
1698   double x0 = (b1 - b2) / (k2 - k1);
1699   //get the radius three ways:
1700   double R1 = sqrt(pow(x1 - x0, 2) + pow(y1 - y0, 2));
1701   double R2 = sqrt(pow(x2 - x0, 2) + pow(y2 - y0, 2));
1702   double R3 = sqrt(pow(x3 - x0, 2) + pow(y3 - y0, 2));
1703   //check if the same
1704   double eps1 = std::abs(R1 / R2 - 1);
1705   double eps2 = std::abs(R3 / R2 - 1);
1706   if (eps1 > 1e-10 || eps2 > 1e-10)
1707     edm::LogVerbatim("Tracklet") << "&&&&&&&&&&&& bad circle! " << R1 << "\t" << R2 << "\t" << R3;
1708 
1709   if (settings_.debugTracklet())
1710     edm::LogVerbatim("Tracklet") << "phimin_: " << phimin_ << " phimax_: " << phimax_;
1711   //results
1712   rinv = 1. / R1;
1713   phi0 = 0.5 * M_PI + atan2(y0, x0);
1714 
1715   phi0 -= phimin_;
1716 
1717   d0 = -R1 + sqrt(x0 * x0 + y0 * y0);
1718   //sign of rinv:
1719   double dphi = reco::reduceRange(phi3 - atan2(y0, x0));
1720   if (dphi < 0) {
1721     rinv = -rinv;
1722     d0 = -d0;
1723     phi0 = phi0 + M_PI;
1724   }
1725   phi0 = angle0to2pi::make0To2pi(phi0);
1726 
1727   //now in RZ:
1728   //turning angle
1729   double beta1 = reco::reduceRange(atan2(y1 - y0, x1 - x0) - atan2(-y0, -x0));
1730   double beta2 = reco::reduceRange(atan2(y2 - y0, x2 - x0) - atan2(-y0, -x0));
1731   double beta3 = reco::reduceRange(atan2(y3 - y0, x3 - x0) - atan2(-y0, -x0));
1732 
1733   double t12 = (z2 - z1) / std::abs(beta2 - beta1) / R1;
1734   double z12 = (z1 * beta2 - z2 * beta1) / (beta2 - beta1);
1735   double t13 = (z3 - z1) / std::abs(beta3 - beta1) / R1;
1736   double z13 = (z1 * beta3 - z3 * beta1) / (beta3 - beta1);
1737 
1738   if (take3 > 0) {
1739     //take 13 (large lever arm)
1740     t = t13;
1741     z0 = z13;
1742   } else {
1743     //take 12 (pixel layers)
1744     t = t12;
1745     z0 = z12;
1746   }
1747 
1748   for (unsigned int i = 0; i < toR_.size(); i++) {
1749     exactproj(toR_[i], rinv, phi0, d0, t, z0, sqrt(x0 * x0 + y0 * y0), phiproj[i], zproj[i], phider[i], zder[i]);
1750   }
1751 
1752   for (unsigned int i = 0; i < toZ_.size(); i++) {
1753     exactprojdisk(toZ_[i], rinv, phi0, d0, t, z0, x0, y0, phiprojdisk[i], rprojdisk[i], phiderdisk[i], rderdisk[i]);
1754   }
1755 
1756   if (settings_.debugTracklet())
1757     edm::LogVerbatim("Tracklet") << "exact tracklet: " << rinv << " " << phi0 << " " << t << " " << z0 << " " << d0;
1758 }
1759 
1760 void TrackletCalculatorDisplaced::approxproj(double halfRinv,
1761                                              double phi0,
1762                                              double d0,
1763                                              double t,
1764                                              double z0,
1765                                              double halfRinv_0,
1766                                              double d0_0,  // zeroeth order result for higher order terms calculation
1767                                              double rmean,
1768                                              double& phiproj,
1769                                              double& phiprojder,
1770                                              double& zproj,
1771                                              double& zprojder) {
1772   if (std::abs(2.0 * halfRinv) > settings_.rinvcut() || std::abs(z0) > settings_.disp_z0cut() ||
1773       std::abs(d0) > settings_.maxd0()) {
1774     phiproj = 0.0;
1775     return;
1776   }
1777   double rmeanInv = 1.0 / rmean;
1778 
1779   phiproj = phi0 + rmean * (-halfRinv + 2.0 * d0_0 * halfRinv_0 * halfRinv_0) +
1780             rmeanInv * (-d0 + halfRinv_0 * d0_0 * d0_0) + sixth * pow(-rmean * halfRinv_0 - rmeanInv * d0_0, 3);
1781   phiprojder = -halfRinv + d0 * rmeanInv * rmeanInv;  //removed all high terms
1782 
1783   zproj = z0 + t * rmean - 0.5 * rmeanInv * t * d0_0 * d0_0 - t * rmean * halfRinv * d0 +
1784           sixth * pow(rmean, 3) * t * halfRinv_0 * halfRinv_0;
1785   zprojder = t;  // removed all high terms
1786 
1787   phiproj = angle0to2pi::make0To2pi(phiproj);
1788 
1789   if (settings_.debugTracklet())
1790     edm::LogVerbatim("Tracklet") << "approx proj layer at " << rmean << " : " << phiproj << " " << zproj << endl;
1791 }
1792 
1793 void TrackletCalculatorDisplaced::approxprojdisk(double halfRinv,
1794                                                  double phi0,
1795                                                  double d0,
1796                                                  double t,
1797                                                  double z0,
1798                                                  double halfRinv_0,
1799                                                  double d0_0,  // zeroeth order result for higher order terms calculation
1800                                                  double zmean,
1801                                                  double& phiproj,
1802                                                  double& phiprojder,
1803                                                  double& rproj,
1804                                                  double& rprojder) {
1805   if (std::abs(2.0 * halfRinv) > settings_.rinvcut() || std::abs(z0) > settings_.disp_z0cut() ||
1806       std::abs(d0) > settings_.maxd0()) {
1807     phiproj = 0.0;
1808     return;
1809   }
1810 
1811   if (t < 0)
1812     zmean = -zmean;
1813 
1814   double zmeanInv = 1.0 / zmean, rstar = (zmean - z0) / t,
1815          epsilon = 0.5 * zmeanInv * zmeanInv * d0_0 * d0_0 * t * t + halfRinv * d0 -
1816                    sixth * rstar * rstar * halfRinv_0 * halfRinv_0;
1817 
1818   rproj = rstar * (1 + epsilon);
1819   rprojder = 1 / t;
1820 
1821   double A = rproj * halfRinv;
1822   double B = -d0 * t * zmeanInv * (1 + z0 * zmeanInv) * (1 - epsilon);
1823   double C = -d0 * halfRinv;
1824   double A_0 = rproj * halfRinv_0;
1825   double B_0 = -d0_0 * t * zmeanInv * (1 + z0 * zmeanInv) * (1 - epsilon);
1826   // double C_0 = -d0_0 * halfRinv_0;
1827 
1828   phiproj = phi0 - A + B * (1 + C - 2 * A_0 * A_0) + sixth * pow(-A_0 + B_0, 3);
1829   phiprojder = -halfRinv / t - d0 * t * t * zmeanInv * zmeanInv;
1830 
1831   phiproj = angle0to2pi::make0To2pi(phiproj);
1832 
1833   if (settings_.debugTracklet())
1834     edm::LogVerbatim("Tracklet") << "approx proj disk at" << zmean << " : " << phiproj << " " << rproj << endl;
1835 }
1836 
1837 void TrackletCalculatorDisplaced::approxtracklet(double r1,
1838                                                  double z1,
1839                                                  double phi1,
1840                                                  double r2,
1841                                                  double z2,
1842                                                  double phi2,
1843                                                  double r3,
1844                                                  double z3,
1845                                                  double phi3,
1846                                                  bool take3,
1847                                                  unsigned ndisks,
1848                                                  double& rinv,
1849                                                  double& phi0,
1850                                                  double& d0,
1851                                                  double& t,
1852                                                  double& z0,
1853                                                  double phiproj[4],
1854                                                  double zproj[4],
1855                                                  double phider[4],
1856                                                  double zder[4],
1857                                                  double phiprojdisk[5],
1858                                                  double rprojdisk[5],
1859                                                  double phiderdisk[5],
1860                                                  double rderdisk[5]) {
1861   double a = 1.0 / ((r1 - r2) * (r1 - r3));
1862   double b = 1.0 / ((r1 - r2) * (r2 - r3));
1863   double c = 1.0 / ((r1 - r3) * (r2 - r3));
1864 
1865   // first iteration in r-phi plane
1866   double halfRinv_0 = -phi1 * r1 * a + phi2 * r2 * b - phi3 * r3 * c;
1867   double d0_0 = r1 * r2 * r3 * (-phi1 * a + phi2 * b - phi3 * c);
1868 
1869   // corrections to phi1, phi2, and phi3
1870   double r = r2, z = z2;
1871   if (take3)
1872     r = r3, z = z3;
1873 
1874   double d0OverR1 = d0_0 * rzmeanInv_[0] * (ndisks > 2 ? std::abs((z - z1) / (r - r1)) : 1.0);
1875   double d0OverR2 = d0_0 * rzmeanInv_[1] * (ndisks > 1 ? std::abs((z - z1) / (r - r1)) : 1.0);
1876   double d0OverR3 = d0_0 * rzmeanInv_[2] * (ndisks > 0 ? std::abs((z - z1) / (r - r1)) : 1.0);
1877 
1878   double d0OverR = d0OverR2;
1879   if (take3)
1880     d0OverR = d0OverR3;
1881 
1882   double c1 = d0_0 * halfRinv_0 * d0OverR1 + 2.0 * d0_0 * halfRinv_0 * r1 * halfRinv_0 +
1883               sixth * pow(-r1 * halfRinv_0 - d0OverR1, 3);
1884   double c2 = d0_0 * halfRinv_0 * d0OverR2 + 2.0 * d0_0 * halfRinv_0 * r2 * halfRinv_0 +
1885               sixth * pow(-r2 * halfRinv_0 - d0OverR2, 3);
1886   double c3 = d0_0 * halfRinv_0 * d0OverR3 + 2.0 * d0_0 * halfRinv_0 * r3 * halfRinv_0 +
1887               sixth * pow(-r3 * halfRinv_0 - d0OverR3, 3);
1888 
1889   double phi1c = phi1 - c1;
1890   double phi2c = phi2 - c2;
1891   double phi3c = phi3 - c3;
1892 
1893   // second iteration in r-phi plane
1894   double halfRinv = -phi1c * r1 * a + phi2c * r2 * b - phi3c * r3 * c;
1895   phi0 = -phi1c * r1 * (r2 + r3) * a + phi2c * r2 * (r1 + r3) * b - phi3c * r3 * (r1 + r2) * c;
1896   d0 = r1 * r2 * r3 * (-phi1c * a + phi2c * b - phi3c * c);
1897 
1898   t = ((z - z1) / (r - r1)) *
1899       (1. + d0 * halfRinv - 0.5 * d0OverR1 * d0OverR - sixth * (r1 * r1 + r2 * r2 + r1 * r2) * halfRinv_0 * halfRinv_0);
1900   z0 = z1 - t * r1 * (1.0 - d0_0 * halfRinv_0 - 0.5 * d0OverR1 * d0OverR1 + sixth * r1 * r1 * halfRinv_0 * halfRinv_0);
1901 
1902   rinv = 2.0 * halfRinv;
1903   phi0 += -phimin_;
1904 
1905   phi0 = angle0to2pi::make0To2pi(phi0);
1906 
1907   for (unsigned int i = 0; i < toR_.size(); i++) {
1908     approxproj(halfRinv,
1909                phi0,
1910                d0,
1911                t,
1912                z0,
1913                halfRinv_0,
1914                d0_0,  // added _0 version for high term calculations
1915                toR_.at(i),
1916                phiproj[i],
1917                phider[i],
1918                zproj[i],
1919                zder[i]);
1920   }
1921 
1922   for (unsigned int i = 0; i < toZ_.size(); i++) {
1923     approxprojdisk(halfRinv,
1924                    phi0,
1925                    d0,
1926                    t,
1927                    z0,
1928                    halfRinv_0,
1929                    d0_0,  // added _0 version for high term calculations
1930                    toZ_.at(i),
1931                    phiprojdisk[i],
1932                    phiderdisk[i],
1933                    rprojdisk[i],
1934                    rderdisk[i]);
1935   }
1936 
1937   if (std::abs(rinv) > settings_.rinvcut() || std::abs(z0) > settings_.disp_z0cut() ||
1938       std::abs(d0) > settings_.maxd0()) {
1939     phi0 = 0.0;
1940     return;
1941   }
1942 
1943   if (settings_.debugTracklet())
1944     edm::LogVerbatim("Tracklet") << "TCD approx tracklet: " << rinv << " " << phi0 << " " << t << " " << z0 << " " << d0
1945                                  << endl;
1946 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!