File indexing completed on 2025-03-23 16:00:17
0001 #include "L1Trigger/TrackFindingTracklet/interface/TrackletLUT.h"
0002 #include "L1Trigger/TrackFindingTracklet/interface/Util.h"
0003 #include "L1Trigger/TrackFindingTracklet/interface/Settings.h"
0004 #include "L1Trigger/TrackFindingTracklet/interface/TrackletConfigBuilder.h"
0005 #include "L1Trigger/TrackTrigger/interface/Setup.h"
0006
0007 #include <filesystem>
0008
0009 using namespace std;
0010 using namespace trklet;
0011
0012 TrackletLUT::TrackletLUT(const Settings& settings)
0013 : settings_(settings), setup_(settings.setup()), nbits_(0), positive_(true) {}
0014
0015 std::vector<const tt::SensorModule*> TrackletLUT::getSensorModules(
0016 unsigned int layerdisk, bool isPS, std::array<double, 2> tan_range, unsigned int nzbins, unsigned int zbin) {
0017
0018
0019
0020
0021
0022
0023 bool use_tan_range = !(tan_range[0] == -1 and tan_range[1] == -1);
0024 bool use_zbins = (nzbins > 1);
0025
0026 bool barrel = layerdisk < N_LAYER;
0027
0028 int layerId = barrel ? layerdisk + 1 : layerdisk + N_LAYER - 1;
0029
0030 std::vector<const tt::SensorModule*> sensorModules;
0031
0032 double z0 = settings_.z0cut();
0033
0034 for (auto& sm : setup_->sensorModules()) {
0035 if (sm.layerId() != layerId || sm.z() < 0 || sm.psModule() != isPS) {
0036 continue;
0037 }
0038
0039 if (use_tan_range) {
0040 const double term = (sm.numColumns() / 2 - 0.5) * sm.pitchCol();
0041 double rmin = sm.r() - term * std::abs(sm.sinTilt());
0042 double rmax = sm.r() + term * std::abs(sm.sinTilt());
0043
0044 double zmin = std::abs(sm.z()) - term * std::abs(sm.cosTilt());
0045 double zmax = std::abs(sm.z()) + term * std::abs(sm.cosTilt());
0046
0047
0048 double mod_tan_max = tan_theta(rmin, zmax, z0, false);
0049 double mod_tan_min = tan_theta(rmax, zmin, z0, true);
0050
0051 if (mod_tan_max >= tan_range[0] && mod_tan_min <= tan_range[1]) {
0052 sensorModules.push_back(&sm);
0053 }
0054 } else if (use_zbins) {
0055 assert(layerdisk < 3);
0056
0057 if (nzbins == 2) {
0058 bool useFlat = (zbin == 0);
0059 bool isFlat = (sm.tilt() == 0);
0060
0061 if (useFlat and isFlat)
0062 sensorModules.push_back(&sm);
0063 else if (!useFlat and !isFlat)
0064 sensorModules.push_back(&sm);
0065 } else if (nzbins == 13) {
0066 if (sm.ringId(setup_) == zbin)
0067 sensorModules.push_back(&sm);
0068 } else {
0069 throw cms::Exception("Unspecified number of z bins");
0070 }
0071 } else {
0072 sensorModules.push_back(&sm);
0073 }
0074 }
0075
0076
0077 static constexpr double delta = 1.e-3;
0078 auto smallerR = [](const tt::SensorModule* lhs, const tt::SensorModule* rhs) { return lhs->r() < rhs->r(); };
0079 auto smallerZ = [](const tt::SensorModule* lhs, const tt::SensorModule* rhs) { return lhs->z() < rhs->z(); };
0080 auto equalRZ = [](const tt::SensorModule* lhs, const tt::SensorModule* rhs) {
0081 return abs(lhs->r() - rhs->r()) < delta && abs(lhs->z() - rhs->z()) < delta;
0082 };
0083 stable_sort(sensorModules.begin(), sensorModules.end(), smallerR);
0084 stable_sort(sensorModules.begin(), sensorModules.end(), smallerZ);
0085 sensorModules.erase(unique(sensorModules.begin(), sensorModules.end(), equalRZ), sensorModules.end());
0086
0087 return sensorModules;
0088 }
0089
0090 std::array<double, 2> TrackletLUT::getTanRange(const std::vector<const tt::SensorModule*>& sensorModules) {
0091
0092
0093 std::array<double, 2> tan_range = {{2147483647, 0}};
0094
0095 double z0 = settings_.z0cut();
0096
0097 for (auto sm : sensorModules) {
0098 const double term = (sm->numColumns() / 2 - 0.5) * sm->pitchCol();
0099 double rmin = sm->r() - term * std::abs(sm->sinTilt());
0100 double rmax = sm->r() + term * std::abs(sm->sinTilt());
0101
0102 double zmin = std::abs(sm->z()) - term * sm->cosTilt();
0103 double zmax = std::abs(sm->z()) + term * sm->cosTilt();
0104
0105 double mod_tan_max = tan_theta(rmin, zmax, z0, true);
0106 double mod_tan_min = tan_theta(rmax, zmin, z0, false);
0107
0108 if (mod_tan_min < tan_range[0])
0109 tan_range[0] = mod_tan_min;
0110 if (mod_tan_max > tan_range[1])
0111 tan_range[1] = mod_tan_max;
0112 }
0113 return tan_range;
0114 }
0115
0116 std::vector<std::array<double, 2>> TrackletLUT::getBendCut(unsigned int layerdisk,
0117 const std::vector<const tt::SensorModule*>& sensorModules,
0118 bool isPS,
0119 double FEbendcut) {
0120
0121
0122
0123
0124
0125 unsigned int bendbits = isPS ? 3 : 4;
0126
0127 std::vector<std::array<double, 2>> bendpars;
0128 std::vector<std::array<double, 2>> bendminmax;
0129
0130
0131 for (int i = 0; i < 1 << bendbits; i++) {
0132 bendpars.push_back({{99, 0}});
0133 bendminmax.push_back({{99, -99}});
0134 }
0135
0136
0137 for (auto sm : sensorModules) {
0138 int window = sm->windowSize();
0139 const vector<double>& encodingBend = setup_->encodingBend(window, isPS);
0140
0141
0142 for (int ibend = 0; ibend <= 2 * window; ibend++) {
0143 int FEbend = ibend - window;
0144 double BEbend = setup_->stubAlgorithm()->degradeBend(isPS, window, FEbend);
0145
0146 const auto pos = std::find(encodingBend.begin(), encodingBend.end(), std::abs(BEbend));
0147 int bend = std::signbit(BEbend) ? (1 << bendbits) - distance(encodingBend.begin(), pos)
0148 : distance(encodingBend.begin(), pos);
0149
0150 double bendmin = FEbend / 2.0 - FEbendcut;
0151 double bendmax = FEbend / 2.0 + FEbendcut;
0152
0153
0154 double z_mod[2];
0155 double r_mod[2];
0156
0157 z_mod[0] = std::abs(sm->z()) + (sm->numColumns() / 2 - 0.5) * sm->pitchCol() * sm->cosTilt();
0158 z_mod[1] = std::abs(sm->z()) - (sm->numColumns() / 2 - 0.5) * sm->pitchCol() * sm->cosTilt();
0159
0160 r_mod[0] = sm->r() - (sm->numColumns() / 2 - 0.5) * sm->pitchCol() * std::abs(sm->sinTilt());
0161 r_mod[1] = sm->r() + (sm->numColumns() / 2 - 0.5) * sm->pitchCol() * std::abs(sm->sinTilt());
0162
0163 for (int i = 0; i < 2; i++) {
0164 double CF = std::abs(sm->sinTilt()) * (z_mod[i] / r_mod[i]) + sm->cosTilt();
0165
0166 double cbendmin =
0167 convertFEBend(bendmin, sm->sep(), settings_.sensorSpacing2S(), CF, (layerdisk < N_LAYER), r_mod[i]);
0168 double cbendmax =
0169 convertFEBend(bendmax, sm->sep(), settings_.sensorSpacing2S(), CF, (layerdisk < N_LAYER), r_mod[i]);
0170
0171 if (cbendmin < bendminmax[bend][0])
0172 bendminmax.at(bend)[0] = cbendmin;
0173 if (cbendmax > bendminmax[bend][1])
0174 bendminmax.at(bend)[1] = cbendmax;
0175 }
0176 }
0177 }
0178
0179 for (int i = 0; i < 1 << bendbits; i++) {
0180 double mid = (bendminmax[i][1] + bendminmax[i][0]) / 2;
0181 double cut = (bendminmax[i][1] - bendminmax[i][0]) / 2;
0182
0183 bendpars[i][0] = mid;
0184 bendpars[i][1] = cut;
0185 }
0186
0187 return bendpars;
0188 }
0189
0190 void TrackletLUT::initmatchcut(unsigned int layerdisk, MatchType type, unsigned int region) {
0191 char cregion = 'A' + region;
0192
0193 for (unsigned int iSeed = 0; iSeed < N_SEED; iSeed++) {
0194 if (type == barrelphi) {
0195 table_.push_back(settings_.rphimatchcut(iSeed, layerdisk) / (settings_.kphi1() * settings_.rmean(layerdisk)));
0196 }
0197 if (type == barrelz) {
0198 table_.push_back(settings_.zmatchcut(iSeed, layerdisk) / settings_.kz());
0199 }
0200 if (type == diskPSphi) {
0201 table_.push_back(settings_.rphicutPS(iSeed, layerdisk - N_LAYER) / (settings_.kphi() * settings_.kr()));
0202 }
0203 if (type == disk2Sphi) {
0204 table_.push_back(settings_.rphicut2S(iSeed, layerdisk - N_LAYER) / (settings_.kphi() * settings_.kr()));
0205 }
0206 if (type == disk2Sr) {
0207 table_.push_back(settings_.rcut2S(iSeed, layerdisk - N_LAYER) / settings_.krprojshiftdisk());
0208 }
0209 if (type == diskPSr) {
0210 table_.push_back(settings_.rcutPS(iSeed, layerdisk - N_LAYER) / settings_.krprojshiftdisk());
0211 }
0212 }
0213 if (type == alphainner) {
0214 for (unsigned int i = 0; i < N_DSS_MOD * 2; i++) {
0215 table_.push_back((1 << settings_.alphashift()) * settings_.krprojshiftdisk() * settings_.half2SmoduleWidth() /
0216 (1 << (settings_.nbitsalpha() - 1)) / (settings_.rDSSinner(i) * settings_.rDSSinner(i)) /
0217 settings_.kphi());
0218 }
0219 }
0220 if (type == alphaouter) {
0221 for (unsigned int i = 0; i < N_DSS_MOD * 2; i++) {
0222 table_.push_back((1 << settings_.alphashift()) * settings_.krprojshiftdisk() * settings_.half2SmoduleWidth() /
0223 (1 << (settings_.nbitsalpha() - 1)) / (settings_.rDSSouter(i) * settings_.rDSSouter(i)) /
0224 settings_.kphi());
0225 }
0226 }
0227 if (type == rSSinner) {
0228 for (unsigned int i = 0; i < N_DSS_MOD * 2; i++) {
0229 table_.push_back(settings_.rDSSinner(i) / settings_.kr());
0230 }
0231 }
0232 if (type == rSSouter) {
0233 for (unsigned int i = 0; i < N_DSS_MOD * 2; i++) {
0234 table_.push_back(settings_.rDSSouter(i) / settings_.kr());
0235 }
0236 }
0237
0238 name_ = settings_.combined() ? "MP_" : "MC_";
0239
0240 if (type == barrelphi) {
0241 nbits_ = 10;
0242 name_ += TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_phicut.tab";
0243 }
0244 if (type == barrelz) {
0245 nbits_ = 9;
0246 name_ += TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_zcut.tab";
0247 }
0248 if (type == diskPSphi) {
0249 nbits_ = 19;
0250 name_ += TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_PSphicut.tab";
0251 }
0252 if (type == disk2Sphi) {
0253 nbits_ = 19;
0254 name_ += TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_2Sphicut.tab";
0255 }
0256 if (type == disk2Sr) {
0257 nbits_ = 7;
0258 name_ += TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_2Srcut.tab";
0259 }
0260 if (type == diskPSr) {
0261 nbits_ = 2;
0262 name_ += TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_PSrcut.tab";
0263 }
0264 if (type == alphainner) {
0265 nbits_ = 8;
0266 name_ += TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_alphainner.tab";
0267 }
0268 if (type == alphaouter) {
0269 nbits_ = 8;
0270 name_ += TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_alphaouter.tab";
0271 }
0272 if (type == rSSinner) {
0273 nbits_ = 15;
0274 name_ += TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_rDSSinner.tab";
0275 }
0276 if (type == rSSouter) {
0277 nbits_ = 15;
0278 name_ += TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_rDSSouter.tab";
0279 }
0280
0281 positive_ = false;
0282
0283 writeTable();
0284 }
0285
0286 void TrackletLUT::initTPlut(bool fillInner,
0287 unsigned int iSeed,
0288 unsigned int layerdisk1,
0289 unsigned int layerdisk2,
0290 unsigned int nbitsfinephidiff,
0291 unsigned int iTP) {
0292
0293 int nfinephibins = settings_.nallstubs(layerdisk2) * settings_.nvmte(1, iSeed) * (1 << settings_.nfinephi(1, iSeed));
0294 double dfinephi = settings_.dphisectorHG() / nfinephibins;
0295
0296 int outerrbits = 3;
0297
0298 if (iSeed == Seed::L1L2 || iSeed == Seed::L2L3 || iSeed == Seed::L3L4 || iSeed == Seed::L5L6) {
0299 outerrbits = 0;
0300 }
0301
0302 int outerrbins = (1 << outerrbits);
0303
0304 double dphi[2];
0305 double router[2];
0306
0307 bool isPSinner;
0308 bool isPSouter;
0309
0310 if (iSeed == Seed::L3L4) {
0311 isPSinner = true;
0312 isPSouter = false;
0313 } else if (iSeed == Seed::L5L6) {
0314 isPSinner = false;
0315 isPSouter = false;
0316 } else {
0317 isPSinner = true;
0318 isPSouter = true;
0319 }
0320
0321 unsigned int nbendbitsinner = isPSinner ? N_BENDBITS_PS : N_BENDBITS_2S;
0322 unsigned int nbendbitsouter = isPSouter ? N_BENDBITS_PS : N_BENDBITS_2S;
0323
0324 double z0 = settings_.z0cut();
0325
0326 int nbinsfinephidiff = (1 << nbitsfinephidiff);
0327
0328 for (int iphibin = 0; iphibin < nbinsfinephidiff; iphibin++) {
0329 int iphidiff = iphibin;
0330 if (iphibin >= nbinsfinephidiff / 2) {
0331 iphidiff = iphibin - nbinsfinephidiff;
0332 }
0333
0334
0335 double deltaphi = 1.5;
0336 dphi[0] = (iphidiff - deltaphi) * dfinephi;
0337 dphi[1] = (iphidiff + deltaphi) * dfinephi;
0338 for (int irouterbin = 0; irouterbin < outerrbins; irouterbin++) {
0339 if (iSeed == Seed::D1D2 || iSeed == Seed::D3D4 || iSeed == Seed::L1D1 || iSeed == Seed::L2D1) {
0340 router[0] =
0341 settings_.rmindiskvm() + irouterbin * (settings_.rmaxdiskvm() - settings_.rmindiskvm()) / outerrbins;
0342 router[1] =
0343 settings_.rmindiskvm() + (irouterbin + 1) * (settings_.rmaxdiskvm() - settings_.rmindiskvm()) / outerrbins;
0344 } else {
0345 router[0] = settings_.rmean(layerdisk2);
0346 router[1] = settings_.rmean(layerdisk2);
0347 }
0348
0349
0350 std::vector<std::array<double, 2>> bend_cuts_inner;
0351 std::vector<std::array<double, 2>> bend_cuts_outer;
0352
0353 if (settings_.useCalcBendCuts) {
0354 std::vector<const tt::SensorModule*> sminner;
0355 std::vector<const tt::SensorModule*> smouter;
0356
0357 if (iSeed == Seed::L1L2 || iSeed == Seed::L2L3 || iSeed == Seed::L3L4 || iSeed == Seed::L5L6) {
0358 double outer_tan_max = tan_theta(settings_.rmean(layerdisk2), settings_.zlength(), z0, true);
0359 std::array<double, 2> tan_range = {{0, outer_tan_max}};
0360
0361 smouter = getSensorModules(layerdisk2, isPSouter, tan_range);
0362 sminner = getSensorModules(layerdisk1, isPSinner, tan_range);
0363
0364 } else if (iSeed == Seed::L1D1 || iSeed == Seed::L2D1) {
0365 double outer_tan_min = tan_theta(router[1], settings_.zmindisk(layerdisk2 - N_LAYER), z0, false);
0366 double outer_tan_max = tan_theta(router[0], settings_.zmaxdisk(layerdisk2 - N_LAYER), z0, true);
0367
0368 smouter = getSensorModules(layerdisk2, isPSouter, {{outer_tan_min, outer_tan_max}});
0369 std::array<double, 2> tan_range = getTanRange(smouter);
0370 sminner = getSensorModules(layerdisk1, isPSinner, tan_range);
0371
0372 } else {
0373
0374 double outer_tan_min = tan_theta(router[1], settings_.zmindisk(layerdisk2 - N_LAYER), z0, false);
0375 double outer_tan_max = tan_theta(router[0], settings_.zmaxdisk(layerdisk2 - N_LAYER), z0, true);
0376
0377 smouter = getSensorModules(layerdisk2, isPSouter, {{outer_tan_min, outer_tan_max}});
0378
0379 std::array<double, 2> tan_range = getTanRange(smouter);
0380 sminner = getSensorModules(layerdisk1, isPSinner, tan_range);
0381 }
0382
0383 bend_cuts_inner = getBendCut(layerdisk1, sminner, isPSinner, settings_.bendcutTE(iSeed, true));
0384 bend_cuts_outer = getBendCut(layerdisk2, smouter, isPSouter, settings_.bendcutTE(iSeed, false));
0385
0386 } else {
0387 for (int ibend = 0; ibend < (1 << nbendbitsinner); ibend++) {
0388 double mid = settings_.benddecode(ibend, layerdisk1, isPSinner);
0389 double cut = settings_.bendcutte(ibend, layerdisk1, isPSinner);
0390 bend_cuts_inner.push_back({{mid, cut}});
0391 }
0392 for (int ibend = 0; ibend < (1 << nbendbitsouter); ibend++) {
0393 double mid = settings_.benddecode(ibend, layerdisk2, isPSouter);
0394 double cut = settings_.bendcutte(ibend, layerdisk2, isPSouter);
0395 bend_cuts_outer.push_back({{mid, cut}});
0396 }
0397 }
0398
0399 double bendinnermin = 20.0;
0400 double bendinnermax = -20.0;
0401 double bendoutermin = 20.0;
0402 double bendoutermax = -20.0;
0403 double rinvmin = 1.0;
0404 double rinvmax = -1.0;
0405 double absrinvmin = 1.0;
0406
0407 for (int i2 = 0; i2 < 2; i2++) {
0408 for (int i3 = 0; i3 < 2; i3++) {
0409 double rinner = 0.0;
0410 if (iSeed == Seed::D1D2 || iSeed == Seed::D3D4) {
0411 rinner = router[i3] * settings_.zmean(layerdisk1 - N_LAYER) / settings_.zmean(layerdisk2 - N_LAYER);
0412 } else {
0413 rinner = settings_.rmean(layerdisk1);
0414 }
0415 if (settings_.useCalcBendCuts) {
0416 if (rinner >= router[i3])
0417 continue;
0418 }
0419 double rinv1 = (rinner < router[i3]) ? rinv(0.0, -dphi[i2], rinner, router[i3]) : 20.0;
0420 double pitchinner = (rinner < settings_.rcrit()) ? settings_.stripPitch(true) : settings_.stripPitch(false);
0421 double pitchouter =
0422 (router[i3] < settings_.rcrit()) ? settings_.stripPitch(true) : settings_.stripPitch(false);
0423 double abendinner = bendstrip(rinner, rinv1, pitchinner, settings_.sensorSpacing2S());
0424 double abendouter = bendstrip(router[i3], rinv1, pitchouter, settings_.sensorSpacing2S());
0425 if (abendinner < bendinnermin)
0426 bendinnermin = abendinner;
0427 if (abendinner > bendinnermax)
0428 bendinnermax = abendinner;
0429 if (abendouter < bendoutermin)
0430 bendoutermin = abendouter;
0431 if (abendouter > bendoutermax)
0432 bendoutermax = abendouter;
0433 if (std::abs(rinv1) < absrinvmin)
0434 absrinvmin = std::abs(rinv1);
0435 if (rinv1 > rinvmax)
0436 rinvmax = rinv1;
0437 if (rinv1 < rinvmin)
0438 rinvmin = rinv1;
0439 }
0440 }
0441
0442 bool passptcut;
0443 double bendfac;
0444 double rinvcutte = settings_.rinvcutte();
0445
0446 if (settings_.useCalcBendCuts) {
0447 double lowrinvcutte =
0448 rinvcutte / 3;
0449 passptcut = rinvmin < rinvcutte and rinvmax > -rinvcutte;
0450 bendfac = (rinvmin < lowrinvcutte and rinvmax > -lowrinvcutte)
0451 ? 1.05
0452 : 1.0;
0453 } else {
0454 passptcut = absrinvmin < rinvcutte;
0455 bendfac = 1.0;
0456 }
0457
0458 if (fillInner) {
0459 for (int ibend = 0; ibend < (1 << nbendbitsinner); ibend++) {
0460 double bendminfac = (isPSinner and (ibend == 2 or ibend == 3)) ? bendfac : 1.0;
0461 double bendmaxfac = (isPSinner and (ibend == 6 or ibend == 5)) ? bendfac : 1.0;
0462
0463 double mid = bend_cuts_inner.at(ibend)[0];
0464 double cut = bend_cuts_inner.at(ibend)[1];
0465
0466 bool passinner = mid + cut * bendmaxfac > bendinnermin && mid - cut * bendminfac < bendinnermax;
0467
0468 table_.push_back(passinner && passptcut);
0469 }
0470 } else {
0471 for (int ibend = 0; ibend < (1 << nbendbitsouter); ibend++) {
0472 double bendminfac = (isPSouter and (ibend == 2 or ibend == 3)) ? bendfac : 1.0;
0473 double bendmaxfac = (isPSouter and (ibend == 6 or ibend == 5)) ? bendfac : 1.0;
0474
0475 double mid = bend_cuts_outer.at(ibend)[0];
0476 double cut = bend_cuts_outer.at(ibend)[1];
0477
0478 bool passouter = mid + cut * bendmaxfac > bendoutermin && mid - cut * bendminfac < bendoutermax;
0479
0480 table_.push_back(passouter && passptcut);
0481 }
0482 }
0483 }
0484 }
0485
0486 positive_ = false;
0487 nbits_ = 1;
0488 char cTP = 'A' + iTP;
0489
0490 name_ = "TP_" + TrackletConfigBuilder::LayerName(layerdisk1) + TrackletConfigBuilder::LayerName(layerdisk2) + cTP;
0491
0492 if (fillInner) {
0493 name_ += "_stubptinnercut.tab";
0494 } else {
0495 name_ += "_stubptoutercut.tab";
0496 }
0497
0498 writeTable();
0499 }
0500
0501 void TrackletLUT::initTPregionlut(unsigned int iSeed,
0502 unsigned int layerdisk1,
0503 unsigned int layerdisk2,
0504 unsigned int iAllStub,
0505 unsigned int nbitsfinephidiff,
0506 unsigned int nbitsfinephi,
0507 const TrackletLUT& tplutinner,
0508 unsigned int iTP) {
0509 int nirbits = 0;
0510 if (iSeed == Seed::D1D2 || iSeed == Seed::D3D4 || iSeed == Seed::L1D1 || iSeed == Seed::L2D1) {
0511 nirbits = 3;
0512 }
0513
0514 unsigned int nbendbitsinner = 3;
0515
0516 if (iSeed == Seed::L5L6) {
0517 nbendbitsinner = 4;
0518 }
0519
0520 for (int innerfinephi = 0; innerfinephi < (1 << nbitsfinephi); innerfinephi++) {
0521 for (int innerbend = 0; innerbend < (1 << nbendbitsinner); innerbend++) {
0522 for (int ir = 0; ir < (1 << nirbits); ir++) {
0523 unsigned int usereg = 0;
0524 for (unsigned int ireg = 0; ireg < settings_.nvmte(1, iSeed); ireg++) {
0525 bool match = false;
0526 for (int ifinephiouter = 0; ifinephiouter < (1 << settings_.nfinephi(1, iSeed)); ifinephiouter++) {
0527 int outerfinephi = iAllStub * (1 << (nbitsfinephi - settings_.nbitsallstubs(layerdisk2))) +
0528 ireg * (1 << settings_.nfinephi(1, iSeed)) + ifinephiouter;
0529 int idphi = outerfinephi - innerfinephi;
0530 bool inrange = (idphi < (1 << (nbitsfinephidiff - 1))) && (idphi >= -(1 << (nbitsfinephidiff - 1)));
0531 if (idphi < 0)
0532 idphi = idphi + (1 << nbitsfinephidiff);
0533 int idphi1 = idphi;
0534 if (iSeed >= 4)
0535 idphi1 = (idphi << 3) + ir;
0536 int ptinnerindexnew = (idphi1 << nbendbitsinner) + innerbend;
0537 match = match || (inrange && tplutinner.lookup(ptinnerindexnew));
0538 }
0539 if (match) {
0540 usereg = usereg | (1 << ireg);
0541 }
0542 }
0543
0544 table_.push_back(usereg);
0545 }
0546 }
0547 }
0548
0549 positive_ = false;
0550 nbits_ = 8;
0551 char cTP = 'A' + iTP;
0552
0553 name_ = "TP_" + TrackletConfigBuilder::LayerName(layerdisk1) + TrackletConfigBuilder::LayerName(layerdisk2) + cTP +
0554 "_usereg.tab";
0555
0556 writeTable();
0557 }
0558
0559 void TrackletLUT::initteptlut(bool fillInner,
0560 bool fillTEMem,
0561 unsigned int iSeed,
0562 unsigned int layerdisk1,
0563 unsigned int layerdisk2,
0564 unsigned int innerphibits,
0565 unsigned int outerphibits,
0566 double innerphimin,
0567 double innerphimax,
0568 double outerphimin,
0569 double outerphimax,
0570 const std::string& innermem,
0571 const std::string& outermem) {
0572 int outerrbits = 0;
0573 if (iSeed == Seed::D1D2 || iSeed == Seed::D3D4 || iSeed == Seed::L1D1 || iSeed == Seed::L2D1) {
0574 outerrbits = 3;
0575 }
0576
0577 int outerrbins = (1 << outerrbits);
0578 int innerphibins = (1 << innerphibits);
0579 int outerphibins = (1 << outerphibits);
0580
0581 double phiinner[2];
0582 double phiouter[2];
0583 double router[2];
0584
0585 bool isPSinner;
0586 bool isPSouter;
0587
0588 if (iSeed == Seed::L3L4) {
0589 isPSinner = true;
0590 isPSouter = false;
0591 } else if (iSeed == Seed::L5L6) {
0592 isPSinner = false;
0593 isPSouter = false;
0594 } else {
0595 isPSinner = true;
0596 isPSouter = true;
0597 }
0598
0599 unsigned int nbendbitsinner = isPSinner ? N_BENDBITS_PS : N_BENDBITS_2S;
0600 unsigned int nbendbitsouter = isPSouter ? N_BENDBITS_PS : N_BENDBITS_2S;
0601
0602 if (fillTEMem) {
0603 if (fillInner) {
0604 table_.resize((1 << nbendbitsinner), false);
0605 } else {
0606 table_.resize((1 << nbendbitsouter), false);
0607 }
0608 }
0609
0610 double z0 = settings_.z0cut();
0611
0612 for (int irouterbin = 0; irouterbin < outerrbins; irouterbin++) {
0613 if (iSeed == Seed::D1D2 || iSeed == Seed::D3D4 || iSeed == Seed::L1D1 || iSeed == Seed::L2D1) {
0614 router[0] = settings_.rmindiskvm() + irouterbin * (settings_.rmaxdiskvm() - settings_.rmindiskvm()) / outerrbins;
0615 router[1] =
0616 settings_.rmindiskvm() + (irouterbin + 1) * (settings_.rmaxdiskvm() - settings_.rmindiskvm()) / outerrbins;
0617 } else {
0618 router[0] = settings_.rmean(layerdisk2);
0619 router[1] = settings_.rmean(layerdisk2);
0620 }
0621
0622
0623 std::vector<std::array<double, 2>> bend_cuts_inner;
0624 std::vector<std::array<double, 2>> bend_cuts_outer;
0625
0626 if (settings_.useCalcBendCuts) {
0627 std::vector<const tt::SensorModule*> sminner;
0628 std::vector<const tt::SensorModule*> smouter;
0629
0630 if (iSeed == Seed::L1L2 || iSeed == Seed::L2L3 || iSeed == Seed::L3L4 || iSeed == Seed::L5L6) {
0631 double outer_tan_max = tan_theta(settings_.rmean(layerdisk2), settings_.zlength(), z0, true);
0632 std::array<double, 2> tan_range = {{0, outer_tan_max}};
0633
0634 smouter = getSensorModules(layerdisk2, isPSouter, tan_range);
0635 sminner = getSensorModules(layerdisk1, isPSinner, tan_range);
0636
0637 } else if (iSeed == Seed::L1D1 || iSeed == Seed::L2D1) {
0638 double outer_tan_min = tan_theta(router[1], settings_.zmindisk(layerdisk2 - N_LAYER), z0, false);
0639 double outer_tan_max = tan_theta(router[0], settings_.zmaxdisk(layerdisk2 - N_LAYER), z0, true);
0640
0641 smouter = getSensorModules(layerdisk2, isPSouter, {{outer_tan_min, outer_tan_max}});
0642 std::array<double, 2> tan_range = getTanRange(smouter);
0643 sminner = getSensorModules(layerdisk1, isPSinner, tan_range);
0644
0645 } else {
0646
0647 double outer_tan_min = tan_theta(router[1], settings_.zmindisk(layerdisk2 - N_LAYER), z0, false);
0648 double outer_tan_max = tan_theta(router[0], settings_.zmaxdisk(layerdisk2 - N_LAYER), z0, true);
0649
0650 smouter = getSensorModules(layerdisk2, isPSouter, {{outer_tan_min, outer_tan_max}});
0651
0652 std::array<double, 2> tan_range = getTanRange(smouter);
0653 sminner = getSensorModules(layerdisk1, isPSinner, tan_range);
0654 }
0655
0656 bend_cuts_inner = getBendCut(layerdisk1, sminner, isPSinner, settings_.bendcutTE(iSeed, true));
0657 bend_cuts_outer = getBendCut(layerdisk2, smouter, isPSouter, settings_.bendcutTE(iSeed, false));
0658
0659 } else {
0660 for (int ibend = 0; ibend < (1 << nbendbitsinner); ibend++) {
0661 double mid = settings_.benddecode(ibend, layerdisk1, nbendbitsinner == 3);
0662 double cut = settings_.bendcutte(ibend, layerdisk1, nbendbitsinner == 3);
0663 bend_cuts_inner.push_back({{mid, cut}});
0664 }
0665 for (int ibend = 0; ibend < (1 << nbendbitsouter); ibend++) {
0666 double mid = settings_.benddecode(ibend, layerdisk2, nbendbitsouter == 3);
0667 double cut = settings_.bendcutte(ibend, layerdisk2, nbendbitsouter == 3);
0668 bend_cuts_outer.push_back({{mid, cut}});
0669 }
0670 }
0671
0672 for (int iphiinnerbin = 0; iphiinnerbin < innerphibins; iphiinnerbin++) {
0673 phiinner[0] = innerphimin + iphiinnerbin * (innerphimax - innerphimin) / innerphibins;
0674 phiinner[1] = innerphimin + (iphiinnerbin + 1) * (innerphimax - innerphimin) / innerphibins;
0675 for (int iphiouterbin = 0; iphiouterbin < outerphibins; iphiouterbin++) {
0676 phiouter[0] = outerphimin + iphiouterbin * (outerphimax - outerphimin) / outerphibins;
0677 phiouter[1] = outerphimin + (iphiouterbin + 1) * (outerphimax - outerphimin) / outerphibins;
0678
0679 double bendinnermin = 20.0;
0680 double bendinnermax = -20.0;
0681 double bendoutermin = 20.0;
0682 double bendoutermax = -20.0;
0683 double rinvmin = 1.0;
0684 double rinvmax = -1.0;
0685 double absrinvmin = 1.0;
0686
0687 for (int i1 = 0; i1 < 2; i1++) {
0688 for (int i2 = 0; i2 < 2; i2++) {
0689 for (int i3 = 0; i3 < 2; i3++) {
0690 double rinner = 0.0;
0691 if (iSeed == Seed::D1D2 || iSeed == Seed::D3D4) {
0692 rinner = router[i3] * settings_.zmean(layerdisk1 - N_LAYER) / settings_.zmean(layerdisk2 - N_LAYER);
0693 } else {
0694 rinner = settings_.rmean(layerdisk1);
0695 }
0696
0697 if (settings_.useCalcBendCuts) {
0698 if (rinner >= router[i3])
0699 continue;
0700 }
0701
0702 double rinv1 = (rinner < router[i3]) ? -rinv(phiinner[i1], phiouter[i2], rinner, router[i3]) : -20.0;
0703 double pitchinner =
0704 (rinner < settings_.rcrit()) ? settings_.stripPitch(true) : settings_.stripPitch(false);
0705 double pitchouter =
0706 (router[i3] < settings_.rcrit()) ? settings_.stripPitch(true) : settings_.stripPitch(false);
0707
0708 double abendinner = bendstrip(rinner, rinv1, pitchinner, settings_.sensorSpacing2S());
0709 double abendouter = bendstrip(router[i3], rinv1, pitchouter, settings_.sensorSpacing2S());
0710
0711 if (abendinner < bendinnermin)
0712 bendinnermin = abendinner;
0713 if (abendinner > bendinnermax)
0714 bendinnermax = abendinner;
0715 if (abendouter < bendoutermin)
0716 bendoutermin = abendouter;
0717 if (abendouter > bendoutermax)
0718 bendoutermax = abendouter;
0719 if (std::abs(rinv1) < absrinvmin)
0720 absrinvmin = std::abs(rinv1);
0721 if (rinv1 > rinvmax)
0722 rinvmax = rinv1;
0723 if (rinv1 < rinvmin)
0724 rinvmin = rinv1;
0725 }
0726 }
0727 }
0728
0729 double lowrinvcutte = 0.002;
0730
0731 bool passptcut;
0732 double bendfac;
0733
0734 if (settings_.useCalcBendCuts) {
0735 passptcut = rinvmin < settings_.rinvcutte() and rinvmax > -settings_.rinvcutte();
0736 bendfac = (rinvmin < lowrinvcutte and rinvmax > -lowrinvcutte) ? 1.05 : 1.0;
0737 } else {
0738 passptcut = absrinvmin < settings_.rinvcutte();
0739 bendfac = 1.0;
0740 }
0741
0742 if (fillInner) {
0743 for (int ibend = 0; ibend < (1 << nbendbitsinner); ibend++) {
0744 double bendminfac = (isPSinner and (ibend == 2 or ibend == 3)) ? bendfac : 1.0;
0745 double bendmaxfac = (isPSinner and (ibend == 6 or ibend == 5)) ? bendfac : 1.0;
0746
0747 double mid = bend_cuts_inner.at(ibend)[0];
0748 double cut = bend_cuts_inner.at(ibend)[1];
0749
0750 bool passinner = mid + cut * bendmaxfac > bendinnermin && mid - cut * bendminfac < bendinnermax;
0751
0752 if (fillTEMem) {
0753 if (passinner)
0754 table_[ibend] = 1;
0755 } else {
0756 table_.push_back(passinner && passptcut);
0757 }
0758 }
0759 } else {
0760 for (int ibend = 0; ibend < (1 << nbendbitsouter); ibend++) {
0761 double bendminfac = (isPSouter and (ibend == 2 or ibend == 3)) ? bendfac : 1.0;
0762 double bendmaxfac = (isPSouter and (ibend == 6 or ibend == 5)) ? bendfac : 1.0;
0763
0764 double mid = bend_cuts_outer.at(ibend)[0];
0765 double cut = bend_cuts_outer.at(ibend)[1];
0766
0767 bool passouter = mid + cut * bendmaxfac > bendoutermin && mid - cut * bendminfac < bendoutermax;
0768
0769 if (fillTEMem) {
0770 if (passouter)
0771 table_[ibend] = 1;
0772 } else {
0773 table_.push_back(passouter && passptcut);
0774 }
0775 }
0776 }
0777 }
0778 }
0779 }
0780
0781 positive_ = false;
0782 nbits_ = 1;
0783
0784 if (fillTEMem) {
0785 if (fillInner) {
0786 name_ = "VMSTE_" + innermem + "_vmbendcut.tab";
0787 } else {
0788 name_ = "VMSTE_" + outermem + "_vmbendcut.tab";
0789 }
0790 } else {
0791 name_ = "TE_" + innermem.substr(0, innermem.size() - 2) + "_" + outermem.substr(0, outermem.size() - 2);
0792 if (fillInner) {
0793 name_ += "_stubptinnercut.tab";
0794 } else {
0795 name_ += "_stubptoutercut.tab";
0796 }
0797 }
0798 writeTable();
0799 }
0800
0801 void TrackletLUT::initProjectionBend(double k_phider,
0802 unsigned int idisk,
0803 unsigned int nrbits,
0804 unsigned int nphiderbits) {
0805 unsigned int nsignbins = 2;
0806 unsigned int nrbins = 1 << (nrbits);
0807 unsigned int nphiderbins = 1 << (nphiderbits);
0808
0809 for (unsigned int isignbin = 0; isignbin < nsignbins; isignbin++) {
0810 for (unsigned int irbin = 0; irbin < nrbins; irbin++) {
0811 int ir = irbin;
0812 if (ir > (1 << (nrbits - 1)))
0813 ir -= (1 << nrbits);
0814 ir = ir << (settings_.nrbitsstub(N_LAYER) - nrbits);
0815 for (unsigned int iphiderbin = 0; iphiderbin < nphiderbins; iphiderbin++) {
0816 int iphider = iphiderbin;
0817 if (iphider > (1 << (nphiderbits - 1)))
0818 iphider -= (1 << nphiderbits);
0819 iphider = iphider << (settings_.nbitsphiprojderL123() - nphiderbits);
0820
0821 double rproj = ir * settings_.krprojshiftdisk();
0822 double phider = iphider * k_phider;
0823 double t = settings_.zmean(idisk) / rproj;
0824
0825 if (isignbin)
0826 t = -t;
0827
0828 double rinv = -phider * (2.0 * t);
0829
0830 double stripPitch = (rproj < settings_.rcrit()) ? settings_.stripPitch(true) : settings_.stripPitch(false);
0831 double bendproj = bendstrip(rproj, rinv, stripPitch, settings_.sensorSpacing2S());
0832
0833 constexpr double maxbend = (1 << NRINVBITS) - 1;
0834
0835 int ibendproj = 2.0 * bendproj + 0.5 * maxbend;
0836 if (ibendproj < 0)
0837 ibendproj = 0;
0838 if (ibendproj > maxbend)
0839 ibendproj = maxbend;
0840
0841 table_.push_back(ibendproj);
0842 }
0843 }
0844 }
0845
0846 positive_ = false;
0847 nbits_ = 5;
0848 name_ = settings_.combined() ? "MP_" : "PR_";
0849 name_ += "ProjectionBend_" + TrackletConfigBuilder::LayerName(N_LAYER + idisk) + ".tab";
0850
0851 writeTable();
0852 }
0853
0854 void TrackletLUT::initProjectionDiskRadius(int nrbits) {
0855
0856
0857
0858
0859 double roffset = 3.0;
0860
0861 for (unsigned int ir = 0; ir < (1u << nrbits); ir++) {
0862 double r = ir * settings_.rmaxdisk() / (1u << nrbits);
0863
0864 int rbin1 =
0865 (1 << N_RZBITS) * (r - roffset - settings_.rmindiskvm()) / (settings_.rmaxdisk() - settings_.rmindiskvm());
0866 int rbin2 =
0867 (1 << N_RZBITS) * (r + roffset - settings_.rmindiskvm()) / (settings_.rmaxdisk() - settings_.rmindiskvm());
0868
0869 if (rbin1 < 0) {
0870 rbin1 = 0;
0871 }
0872 rbin2 = clamp(rbin2, 0, ((1 << N_RZBITS) - 1));
0873
0874 assert(rbin1 <= rbin2);
0875 assert(rbin2 - rbin1 <= 1);
0876
0877 int d = rbin1 != rbin2;
0878
0879 int finer =
0880 (1 << (N_RZBITS + NFINERZBITS)) *
0881 ((r - settings_.rmindiskvm()) - rbin1 * (settings_.rmaxdisk() - settings_.rmindiskvm()) / (1 << N_RZBITS)) /
0882 (settings_.rmaxdisk() - settings_.rmindiskvm());
0883
0884 finer = clamp(finer, 0, ((1 << (NFINERZBITS + 1)) - 1));
0885
0886
0887 int N_DIFF_FLAG = 1;
0888
0889 int word = (finer << (N_RZBITS + N_DIFF_FLAG)) + (rbin1 << N_DIFF_FLAG) + d;
0890
0891 table_.push_back(word);
0892 }
0893
0894
0895 nbits_ = NFINERZBITS + 1 + N_RZBITS + 1;
0896 positive_ = true;
0897 name_ = "ProjectionDiskRadius.tab";
0898 writeTable();
0899 }
0900
0901 void TrackletLUT::initBendMatch(unsigned int layerdisk) {
0902 unsigned int nrinv = NRINVBITS;
0903 double rinvhalf = 0.5 * ((1 << nrinv) - 1);
0904
0905 bool barrel = layerdisk < N_LAYER;
0906
0907 if (barrel) {
0908 bool isPSmodule = layerdisk < N_PSLAYER;
0909 double stripPitch = settings_.stripPitch(isPSmodule);
0910 unsigned int nbits = isPSmodule ? N_BENDBITS_PS : N_BENDBITS_2S;
0911
0912 std::vector<std::array<double, 2>> bend_cuts;
0913
0914 if (settings_.useCalcBendCuts) {
0915 double bendcutFE = settings_.bendcutME(layerdisk, isPSmodule);
0916 std::vector<const tt::SensorModule*> sm = getSensorModules(layerdisk, isPSmodule);
0917 bend_cuts = getBendCut(layerdisk, sm, isPSmodule, bendcutFE);
0918
0919 } else {
0920 for (unsigned int ibend = 0; ibend < (1u << nbits); ibend++) {
0921 double mid = settings_.benddecode(ibend, layerdisk, isPSmodule);
0922 double cut = settings_.bendcutte(ibend, layerdisk, isPSmodule);
0923 bend_cuts.push_back({{mid, cut}});
0924 }
0925 }
0926
0927 for (unsigned int irinv = 0; irinv < (1u << nrinv); irinv++) {
0928 double rinv = (irinv - rinvhalf) * (1 << (settings_.nbitsrinv() - nrinv)) * settings_.krinvpars();
0929
0930 double projbend = bendstrip(settings_.rmean(layerdisk), rinv, stripPitch, settings_.sensorSpacing2S());
0931 for (unsigned int ibend = 0; ibend < (1u << nbits); ibend++) {
0932 double mid = bend_cuts[ibend][0];
0933 double cut = bend_cuts[ibend][1];
0934
0935 double pass = mid + cut > projbend && mid - cut < projbend;
0936
0937 table_.push_back(pass);
0938 }
0939 }
0940 } else {
0941 std::vector<std::array<double, 2>> bend_cuts_2S;
0942 std::vector<std::array<double, 2>> bend_cuts_PS;
0943
0944 if (settings_.useCalcBendCuts) {
0945 double bendcutFE2S = settings_.bendcutME(layerdisk, false);
0946 std::vector<const tt::SensorModule*> sm2S = getSensorModules(layerdisk, false);
0947 bend_cuts_2S = getBendCut(layerdisk, sm2S, false, bendcutFE2S);
0948
0949 double bendcutFEPS = settings_.bendcutME(layerdisk, true);
0950 std::vector<const tt::SensorModule*> smPS = getSensorModules(layerdisk, true);
0951 bend_cuts_PS = getBendCut(layerdisk, smPS, true, bendcutFEPS);
0952
0953 } else {
0954 for (unsigned int ibend = 0; ibend < (1 << N_BENDBITS_2S); ibend++) {
0955 double mid = settings_.benddecode(ibend, layerdisk, false);
0956 double cut = settings_.bendcutme(ibend, layerdisk, false);
0957 bend_cuts_2S.push_back({{mid, cut}});
0958 }
0959 for (unsigned int ibend = 0; ibend < (1 << N_BENDBITS_PS); ibend++) {
0960 double mid = settings_.benddecode(ibend, layerdisk, true);
0961 double cut = settings_.bendcutme(ibend, layerdisk, true);
0962 bend_cuts_PS.push_back({{mid, cut}});
0963 }
0964 }
0965
0966 for (unsigned int iprojbend = 0; iprojbend < (1u << nrinv); iprojbend++) {
0967 double projbend = 0.5 * (iprojbend - rinvhalf);
0968 for (unsigned int ibend = 0; ibend < (1 << N_BENDBITS_2S); ibend++) {
0969 double mid = bend_cuts_2S[ibend][0];
0970 double cut = bend_cuts_2S[ibend][1];
0971
0972 double pass = mid + cut > projbend && mid - cut < projbend;
0973
0974 table_.push_back(pass);
0975 }
0976 }
0977 for (unsigned int iprojbend = 0; iprojbend < (1u << nrinv); iprojbend++) {
0978 double projbend = 0.5 * (iprojbend - rinvhalf);
0979 for (unsigned int ibend = 0; ibend < (1 << N_BENDBITS_PS); ibend++) {
0980 double mid = bend_cuts_PS[ibend][0];
0981 double cut = bend_cuts_PS[ibend][1];
0982
0983 double pass = mid + cut > projbend && mid - cut < projbend;
0984
0985 table_.push_back(pass);
0986 }
0987 }
0988 }
0989
0990 positive_ = false;
0991 nbits_ = 1;
0992
0993 name_ = "METable_" + TrackletConfigBuilder::LayerName(layerdisk) + ".tab";
0994
0995 writeTable();
0996 }
0997
0998 void TrackletLUT::initVMRTable(unsigned int layerdisk, VMRTableType type, int region, bool combined) {
0999 unsigned int zbits = settings_.vmrlutzbits(layerdisk);
1000 unsigned int rbits = settings_.vmrlutrbits(layerdisk);
1001
1002 unsigned int rbins = (1 << rbits);
1003 unsigned int zbins = (1 << zbits);
1004
1005 double zmin, zmax, rmin, rmax;
1006
1007 if (layerdisk < N_LAYER) {
1008 zmin = -settings_.zlength();
1009 zmax = settings_.zlength();
1010 rmin = settings_.rmean(layerdisk) - settings_.drmax();
1011 rmax = settings_.rmean(layerdisk) + settings_.drmax();
1012 } else {
1013 rmin = 0;
1014 rmax = settings_.rmaxdisk();
1015 zmin = settings_.zmean(layerdisk - N_LAYER) - settings_.dzmax();
1016 zmax = settings_.zmean(layerdisk - N_LAYER) + settings_.dzmax();
1017 }
1018
1019 double dr = (rmax - rmin) / rbins;
1020 double dz = (zmax - zmin) / zbins;
1021
1022 int NBINS = settings_.NLONGVMBINS() * settings_.NLONGVMBINS();
1023
1024 for (unsigned int izbin = 0; izbin < zbins; izbin++) {
1025 for (unsigned int irbin = 0; irbin < rbins; irbin++) {
1026 double r = rmin + (irbin + 0.5) * dr;
1027 double z = zmin + (izbin + 0.5) * dz;
1028
1029
1030
1031
1032
1033
1034 if (settings_.combined() && combined) {
1035 int iznew = izbin - (1 << (zbits - 1));
1036 if (iznew < 0)
1037 iznew += (1 << zbits);
1038 assert(iznew >= 0);
1039 assert(iznew < (1 << zbits));
1040 z = zmin + (iznew + 0.5) * dz;
1041 if (layerdisk < N_LAYER) {
1042 int irnew = irbin - (1 << (rbits - 1));
1043 if (irnew < 0)
1044 irnew += (1 << rbits);
1045 assert(irnew >= 0);
1046 assert(irnew < (1 << rbits));
1047 r = rmin + (irnew + 0.5) * dr;
1048 }
1049 }
1050
1051 unsigned int NRING =
1052 5;
1053 if (layerdisk >= N_LAYER && irbin < 2 * NRING)
1054 r = (layerdisk < N_LAYER + 2) ? settings_.rDSSinner(irbin) : settings_.rDSSouter(irbin);
1055
1056 int bin;
1057 if (layerdisk < N_LAYER) {
1058 double zproj = z * settings_.rmean(layerdisk) / r;
1059 bin = NBINS * (zproj + settings_.zlength()) / (2 * settings_.zlength());
1060 } else {
1061 double rproj = r * settings_.zmean(layerdisk - N_LAYER) / z;
1062 bin = NBINS * (rproj - settings_.rmindiskvm()) / (settings_.rmaxdisk() - settings_.rmindiskvm());
1063 }
1064 if (bin < 0)
1065 bin = 0;
1066 if (bin >= NBINS)
1067 bin = NBINS - 1;
1068
1069 if (type == VMRTableType::me) {
1070 table_.push_back(bin);
1071 }
1072
1073 if (type == VMRTableType::disk) {
1074 if (layerdisk >= N_LAYER) {
1075 double rproj = r * settings_.zmean(layerdisk - N_LAYER) / z;
1076 bin = 0.5 * NBINS * (rproj - settings_.rmindiskvm()) / (settings_.rmaxdiskvm() - settings_.rmindiskvm());
1077
1078 if (bin < 0)
1079 bin = 0;
1080 if (bin >= NBINS / 2)
1081 bin = 0;
1082 table_.push_back(bin);
1083 }
1084 }
1085
1086 if (type == VMRTableType::inner) {
1087 if (layerdisk == LayerDisk::L1 || layerdisk == LayerDisk::L3 || layerdisk == LayerDisk::L5 ||
1088 layerdisk == LayerDisk::D1 || layerdisk == LayerDisk::D3) {
1089 table_.push_back(getVMRLookup(layerdisk + 1, z, r, dz, dr));
1090 }
1091 if (layerdisk == LayerDisk::L2) {
1092 table_.push_back(getVMRLookup(layerdisk + 1, z, r, dz, dr, Seed::L2L3));
1093 }
1094 }
1095
1096 if (type == VMRTableType::inneroverlap) {
1097 if (layerdisk == LayerDisk::L1 || layerdisk == LayerDisk::L2) {
1098 table_.push_back(getVMRLookup(6, z, r, dz, dr, layerdisk + 6));
1099 }
1100 }
1101
1102 if (type == VMRTableType::innerthird) {
1103 if (layerdisk == LayerDisk::L2) {
1104 table_.push_back(getVMRLookup(LayerDisk::D1, z, r, dz, dr, Seed::L2L3D1));
1105 }
1106
1107 if (layerdisk == LayerDisk::L5) {
1108 table_.push_back(getVMRLookup(LayerDisk::L4, z, r, dz, dr));
1109 }
1110
1111 if (layerdisk == LayerDisk::L3) {
1112 table_.push_back(getVMRLookup(LayerDisk::L2, z, r, dz, dr));
1113 }
1114
1115 if (layerdisk == LayerDisk::D1) {
1116 table_.push_back(getVMRLookup(LayerDisk::L2, z, r, dz, dr));
1117 }
1118 }
1119 }
1120 }
1121
1122
1123
1124
1125
1126
1127 if (settings_.combined() && combined) {
1128 if (type == VMRTableType::me) {
1129 nbits_ = 2 * settings_.NLONGVMBITS();
1130 positive_ = false;
1131 name_ = "VMRME_" + TrackletConfigBuilder::LayerName(layerdisk) + ".tab";
1132 }
1133 if (type == VMRTableType::disk) {
1134 nbits_ = 2 * settings_.NLONGVMBITS();
1135 positive_ = false;
1136 name_ = "VMRTE_" + TrackletConfigBuilder::LayerName(layerdisk) + ".tab";
1137 }
1138 if (type == VMRTableType::inner) {
1139 positive_ = true;
1140 nbits_ = 10;
1141 name_ = "TP_" + TrackletConfigBuilder::LayerName(layerdisk) + TrackletConfigBuilder::LayerName(layerdisk + 1) +
1142 ".tab";
1143 }
1144
1145 if (type == VMRTableType::inneroverlap) {
1146 positive_ = true;
1147 nbits_ = 10;
1148 name_ = "TP_" + TrackletConfigBuilder::LayerName(layerdisk) + TrackletConfigBuilder::LayerName(N_LAYER) + ".tab";
1149 }
1150
1151 } else {
1152 if (type == VMRTableType::me) {
1153
1154 if (layerdisk == LayerDisk::L2 || layerdisk == LayerDisk::L3 || layerdisk == LayerDisk::L4 ||
1155 layerdisk == LayerDisk::L6) {
1156 nbits_ = 6;
1157 positive_ = false;
1158 name_ = "VMTableOuter" + TrackletConfigBuilder::LayerName(layerdisk) + ".tab";
1159 writeTable();
1160 }
1161
1162 assert(region >= 0);
1163 char cregion = 'A' + region;
1164 name_ = "VMR_" + TrackletConfigBuilder::LayerName(layerdisk) + "PHI" + cregion + "_finebin.tab";
1165 nbits_ = 6;
1166 positive_ = false;
1167 }
1168
1169 if (type == VMRTableType::inner) {
1170 positive_ = false;
1171 nbits_ = 10;
1172 name_ = "VMTableInner" + TrackletConfigBuilder::LayerName(layerdisk) +
1173 TrackletConfigBuilder::LayerName(layerdisk + 1) + ".tab";
1174 }
1175
1176 if (type == VMRTableType::inneroverlap) {
1177 positive_ = false;
1178 nbits_ = 10;
1179 name_ = "VMTableInner" + TrackletConfigBuilder::LayerName(layerdisk) + TrackletConfigBuilder::LayerName(N_LAYER) +
1180 ".tab";
1181 }
1182
1183 if (type == VMRTableType::disk) {
1184 positive_ = false;
1185 nbits_ = 10;
1186 name_ = "VMTableOuter" + TrackletConfigBuilder::LayerName(layerdisk) + ".tab";
1187 }
1188 }
1189
1190 writeTable();
1191 }
1192
1193 int TrackletLUT::getVMRLookup(unsigned int layerdisk, double z, double r, double dz, double dr, int iseed) const {
1194 double z0cut = settings_.z0cut();
1195
1196 if (layerdisk < N_LAYER) {
1197 double constexpr zcutL2L3 = 52.0;
1198 if (iseed == Seed::L2L3 && std::abs(z) < zcutL2L3)
1199 return -1;
1200
1201 double rmean = settings_.rmean(layerdisk);
1202
1203 double rratio1 = rmean / (r + 0.5 * dr);
1204 double rratio2 = rmean / (r - 0.5 * dr);
1205
1206 double z1 = (z - 0.5 * dz) * rratio1 + z0cut * (rratio1 - 1.0);
1207 double z2 = (z + 0.5 * dz) * rratio1 + z0cut * (rratio1 - 1.0);
1208 double z3 = (z - 0.5 * dz) * rratio2 + z0cut * (rratio2 - 1.0);
1209 double z4 = (z + 0.5 * dz) * rratio2 + z0cut * (rratio2 - 1.0);
1210 double z5 = (z - 0.5 * dz) * rratio1 - z0cut * (rratio1 - 1.0);
1211 double z6 = (z + 0.5 * dz) * rratio1 - z0cut * (rratio1 - 1.0);
1212 double z7 = (z - 0.5 * dz) * rratio2 - z0cut * (rratio2 - 1.0);
1213 double z8 = (z + 0.5 * dz) * rratio2 - z0cut * (rratio2 - 1.0);
1214
1215 double zmin = std::min({z1, z2, z3, z4, z5, z6, z7, z8});
1216 double zmax = std::max({z1, z2, z3, z4, z5, z6, z7, z8});
1217
1218 int NBINS = settings_.NLONGVMBINS() * settings_.NLONGVMBINS();
1219
1220 int zbin1 = NBINS * (zmin + settings_.zlength()) / (2 * settings_.zlength());
1221 int zbin2 = NBINS * (zmax + settings_.zlength()) / (2 * settings_.zlength());
1222
1223 if (zbin1 >= NBINS)
1224 return -1;
1225 if (zbin2 < 0)
1226 return -1;
1227
1228 if (zbin2 >= NBINS)
1229 zbin2 = NBINS - 1;
1230 if (zbin1 < 0)
1231 zbin1 = 0;
1232
1233
1234
1235
1236
1237
1238
1239
1240
1241
1242
1243 int value = zbin1 / 8;
1244 value *= 2;
1245 if (zbin2 / 8 - zbin1 / 8 > 0)
1246 value += 1;
1247 value *= 8;
1248 value += (zbin1 & 7);
1249 assert(value / 8 < 15);
1250 int deltaz = zbin2 - zbin1;
1251 if (deltaz > 7) {
1252 deltaz = 7;
1253 }
1254 assert(deltaz < 8);
1255 value += (deltaz << 7);
1256
1257 return value;
1258
1259 } else {
1260 if (std::abs(z) < 2.0 * z0cut)
1261 return -1;
1262
1263 double zmean = settings_.zmean(layerdisk - N_LAYER);
1264 if (z < 0.0)
1265 zmean = -zmean;
1266
1267 double r1 = (r + 0.5 * dr) * (zmean + z0cut) / (z + 0.5 * dz + z0cut);
1268 double r2 = (r - 0.5 * dr) * (zmean - z0cut) / (z + 0.5 * dz - z0cut);
1269 double r3 = (r + 0.5 * dr) * (zmean + z0cut) / (z - 0.5 * dz + z0cut);
1270 double r4 = (r - 0.5 * dr) * (zmean - z0cut) / (z - 0.5 * dz - z0cut);
1271 double r5 = (r + 0.5 * dr) * (zmean - z0cut) / (z + 0.5 * dz - z0cut);
1272 double r6 = (r - 0.5 * dr) * (zmean + z0cut) / (z + 0.5 * dz + z0cut);
1273 double r7 = (r + 0.5 * dr) * (zmean - z0cut) / (z - 0.5 * dz - z0cut);
1274 double r8 = (r - 0.5 * dr) * (zmean + z0cut) / (z - 0.5 * dz + z0cut);
1275
1276 double rmin = std::min({r1, r2, r3, r4, r5, r6, r7, r8});
1277 double rmax = std::max({r1, r2, r3, r4, r5, r6, r7, r8});
1278
1279 int NBINS = settings_.NLONGVMBINS() * settings_.NLONGVMBINS() / 2;
1280
1281 double rmindisk = settings_.rmindiskvm();
1282 double rmaxdisk = settings_.rmaxdiskvm();
1283
1284 if (iseed == Seed::L1D1)
1285 rmaxdisk = settings_.rmaxdiskl1overlapvm();
1286 if (iseed == Seed::L2D1)
1287 rmindisk = settings_.rmindiskl2overlapvm();
1288 if (iseed == Seed::L2L3D1)
1289 rmaxdisk = settings_.rmaxdisk();
1290
1291 if (rmin > rmaxdisk)
1292 return -1;
1293 if (rmax > rmaxdisk)
1294 rmax = rmaxdisk;
1295
1296 if (rmax < rmindisk)
1297 return -1;
1298 if (rmin < rmindisk)
1299 rmin = rmindisk;
1300
1301 int rbin1 = NBINS * (rmin - settings_.rmindiskvm()) / (settings_.rmaxdiskvm() - settings_.rmindiskvm());
1302 int rbin2 = NBINS * (rmax - settings_.rmindiskvm()) / (settings_.rmaxdiskvm() - settings_.rmindiskvm());
1303
1304 if (iseed == Seed::L2L3D1) {
1305 constexpr double rminspec = 40.0;
1306 rbin1 = NBINS * (rmin - rminspec) / (settings_.rmaxdisk() - rminspec);
1307 rbin2 = NBINS * (rmax - rminspec) / (settings_.rmaxdisk() - rminspec);
1308 }
1309
1310 if (rbin2 >= NBINS)
1311 rbin2 = NBINS - 1;
1312 if (rbin1 < 0)
1313 rbin1 = 0;
1314
1315
1316
1317
1318
1319
1320
1321
1322
1323
1324
1325 bool overlap = iseed == Seed::L1D1 || iseed == Seed::L2D1 || iseed == Seed::L2L3D1;
1326
1327 int value = rbin1 / 8;
1328 if (overlap) {
1329 if (z < 0.0)
1330 value += 4;
1331 }
1332 value *= 2;
1333 if (rbin2 / 8 - rbin1 / 8 > 0)
1334 value += 1;
1335 value *= 8;
1336 value += (rbin1 & 7);
1337 assert(value / 8 < 15);
1338 int deltar = rbin2 - rbin1;
1339 if (deltar > 7)
1340 deltar = 7;
1341 if (overlap) {
1342 value += (deltar << 7);
1343 } else {
1344 value += (deltar << 6);
1345 }
1346
1347 return value;
1348 }
1349 }
1350
1351 void TrackletLUT::initPhiCorrTable(unsigned int layerdisk, unsigned int rbits) {
1352 bool psmodule = layerdisk < N_PSLAYER;
1353
1354 unsigned int bendbits = psmodule ? N_BENDBITS_PS : N_BENDBITS_2S;
1355
1356 unsigned int rbins = (1 << rbits);
1357
1358 double rmean = settings_.rmean(layerdisk);
1359 double drmax = settings_.drmax();
1360
1361 double dr = 2.0 * drmax / rbins;
1362
1363 std::vector<std::array<double, 2>> bend_vals;
1364
1365 if (settings_.useCalcBendCuts) {
1366 std::vector<const tt::SensorModule*> sm = getSensorModules(layerdisk, psmodule);
1367 bend_vals = getBendCut(layerdisk, sm, psmodule);
1368
1369 } else {
1370 for (int ibend = 0; ibend < 1 << bendbits; ibend++) {
1371 bend_vals.push_back({{settings_.benddecode(ibend, layerdisk, layerdisk < N_PSLAYER), 0}});
1372 }
1373 }
1374
1375 for (int ibend = 0; ibend < 1 << bendbits; ibend++) {
1376 for (unsigned int irbin = 0; irbin < rbins; irbin++) {
1377 double bend = -bend_vals[ibend][0];
1378 int value = getphiCorrValue(layerdisk, bend, irbin, rmean, dr, drmax);
1379 table_.push_back(value);
1380 }
1381 }
1382
1383 name_ = "VMPhiCorrL" + std::to_string(layerdisk + 1) + ".tab";
1384 nbits_ = 14;
1385 positive_ = false;
1386
1387 writeTable();
1388 }
1389
1390 int TrackletLUT::getphiCorrValue(
1391 unsigned int layerdisk, double bend, unsigned int irbin, double rmean, double dr, double drmax) const {
1392 bool psmodule = layerdisk < N_PSLAYER;
1393
1394
1395 double Delta = (irbin + 0.5) * dr - drmax;
1396
1397
1398 double drnom = 0.18;
1399 double dphi = (Delta / drnom) * bend * settings_.stripPitch(psmodule) / rmean;
1400
1401 double kphi = psmodule ? settings_.kphi() : settings_.kphi1();
1402
1403 int idphi = dphi / kphi;
1404
1405 return idphi;
1406 }
1407
1408
1409 void TrackletLUT::writeTable() const {
1410 if (name_.empty()) {
1411 return;
1412 }
1413
1414 if (nbits_ == 0) {
1415 throw cms::Exception("LogicError") << "Error in " << __FILE__ << " nbits_ == 0 ";
1416 }
1417
1418 if (!settings_.writeTable()) {
1419 return;
1420 }
1421
1422 ofstream out = openfile(settings_.tablePath(), name_, __FILE__, __LINE__);
1423
1424 out << "{" << endl;
1425 for (unsigned int i = 0; i < table_.size(); i++) {
1426 if (i != 0) {
1427 out << "," << endl;
1428 }
1429
1430 int itable = table_[i];
1431 if (positive_) {
1432 if (table_[i] < 0) {
1433 itable = (1 << nbits_) - 1;
1434 }
1435 }
1436
1437 out << itable;
1438 }
1439 out << endl << "};" << endl;
1440 out.close();
1441
1442 string name = name_;
1443
1444 name[name_.size() - 3] = 'd';
1445 name[name_.size() - 2] = 'a';
1446 name[name_.size() - 1] = 't';
1447
1448 out = openfile(settings_.tablePath(), name, __FILE__, __LINE__);
1449
1450 int width = (nbits_ + 3) / 4;
1451
1452 for (unsigned int i = 0; i < table_.size(); i++) {
1453 int itable = table_[i];
1454 if (positive_) {
1455 if (table_[i] < 0) {
1456 itable = (1 << nbits_) - 1;
1457 }
1458 }
1459
1460 out << uppercase << setfill('0') << setw(width) << hex << itable << dec << endl;
1461 }
1462
1463 out.close();
1464 }
1465
1466 int TrackletLUT::lookup(unsigned int index) const {
1467 if (index >= table_.size()) {
1468 throw cms::Exception("LogicError") << "Error in " << __FILE__ << " index >= size " << index << " " << table_.size()
1469 << " in " << name_;
1470 }
1471 assert(index < table_.size());
1472 return table_[index];
1473 }