File indexing completed on 2022-10-14 01:44:10
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);
0032 if (name1[3] == 'L')
0033 layer_ = name1[4] - '0';
0034 if (name1[3] == 'D')
0035 disk_ = name1[4] - '0';
0036
0037
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
0161 auto* tmp = dynamic_cast<TrackletProjectionsMemory*>(memory);
0162 assert(tmp != nullptr);
0163
0164 unsigned int layerdisk = output[8] - '1';
0165 unsigned int phiregion = output[12] - 'A';
0166
0167 if (output[7] == 'L') {
0168 assert(layerdisk < N_LAYER);
0169 assert(phiregion < trackletprojlayers_[layerdisk].size());
0170
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
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
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
0261 if (innerFPGAStub->layerdisk() < N_LAYER && middleFPGAStub->layerdisk() >= N_LAYER &&
0262 outerFPGAStub->layerdisk() >= N_LAYER) {
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) {
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;
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
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
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
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
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) + (d0approx / settings_.rcrit()));
0595 int phicrit = iphi0 - 2 * irinv - 2 * id0;
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
0626 if (izproj[i] < -(1 << (settings_.nzbitsstub(0) - 1)))
0627 continue;
0628 if (izproj[i] >= (1 << (settings_.nzbitsstub(0) - 1)))
0629 continue;
0630
0631
0632 if (iphiproj[i] >= (1 << settings_.nphibitsstub(N_LAYER - 1)) - 1)
0633 continue;
0634 if (iphiproj[i] <= 0)
0635 continue;
0636
0637
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;
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
0688 if (iphiprojdisk[i] <= 0)
0689 continue;
0690 if (iphiprojdisk[i] >= (1 << settings_.nphibitsstub(0)) - 1)
0691 continue;
0692
0693
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;
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
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
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
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
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) + (d0approx / settings_.rcrit()));
1010 int phicrit = iphi0 - 2 * irinv - 2 * id0;
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
1041 if (izproj[i] < -(1 << (settings_.nzbitsstub(0) - 1)))
1042 continue;
1043 if (izproj[i] >= (1 << (settings_.nzbitsstub(0) - 1)))
1044 continue;
1045
1046
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;
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
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
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
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
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) + (d0approx / settings_.rcrit()));
1406 int phicrit = iphi0 - 2 * irinv - 2 * id0;
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
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
1491 if (iphiprojdisk[i] <= 0)
1492 continue;
1493 if (iphiprojdisk[i] >= (1 << settings_.nphibitsstub(0)) - 1)
1494 continue;
1495
1496
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
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,
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
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
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
1686
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
1697 double y0 = (b1 * k2 - b2 * k1) / (k2 - k1);
1698 double x0 = (b1 - b2) / (k2 - k1);
1699
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
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
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
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
1728
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
1740 t = t13;
1741 z0 = z13;
1742 } else {
1743
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,
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;
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;
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,
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
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 * 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
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
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
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,
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,
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 }