File indexing completed on 2024-09-05 05:04:36
0001 #ifndef L1TRIGGER_L1CALOTRIGGER_PHASE2L1CALOJETEMULATOR_H
0002 #define L1TRIGGER_L1CALOTRIGGER_PHASE2L1CALOJETEMULATOR_H
0003
0004 #include "DataFormats/L1TCalorimeterPhase2/interface/Phase2L1CaloJet.h"
0005 #include <cstdlib>
0006
0007 static constexpr int nBarrelEta = 34;
0008 static constexpr int nBarrelPhi = 72;
0009 static constexpr int nHgcalEta = 36;
0010 static constexpr int nHgcalPhi = 72;
0011 static constexpr int nHfEta = 24;
0012 static constexpr int nHfPhi = 72;
0013 static constexpr int nSTEta = 8;
0014 static constexpr int nSTPhi = 24;
0015 static constexpr int nJets = 10;
0016
0017 namespace gctobj {
0018
0019 class towerMax {
0020 public:
0021 float energy;
0022 int phi;
0023 int eta;
0024 float energyMax;
0025 int phiMax;
0026 int etaMax;
0027 int phiCenter;
0028 int etaCenter;
0029
0030 towerMax() {
0031 energy = 0;
0032 phi = 0;
0033 eta = 0;
0034 energyMax = 0;
0035 phiMax = 0;
0036 etaMax = 0;
0037 phiCenter = 0;
0038 etaCenter = 0;
0039 }
0040 };
0041
0042 class jetInfo {
0043 public:
0044 float seedEnergy;
0045 float energy;
0046 float tauEt;
0047 int phi;
0048 int eta;
0049 float energyMax;
0050 int phiMax;
0051 int etaMax;
0052 int phiCenter;
0053 int etaCenter;
0054
0055 jetInfo() {
0056 seedEnergy = 0;
0057 energy = 0;
0058 tauEt = 0;
0059 phi = 0;
0060 eta = 0;
0061 energyMax = 0;
0062 phiMax = 0;
0063 etaMax = 0;
0064 phiCenter = 0;
0065 etaCenter = 0;
0066 }
0067 };
0068
0069 typedef struct {
0070 float et;
0071 int eta;
0072 int phi;
0073 float towerEt;
0074 int towerEta;
0075 int towerPhi;
0076 int centerEta;
0077 int centerPhi;
0078 } GCTsupertower_t;
0079
0080 typedef struct {
0081 GCTsupertower_t cr[nSTPhi];
0082 } etaStrip_t;
0083
0084 typedef struct {
0085 GCTsupertower_t etaStrip[nSTEta];
0086 } hgcalRegion_t;
0087
0088 typedef struct {
0089 GCTsupertower_t pk[nSTEta];
0090 } etaStripPeak_t;
0091
0092 typedef struct {
0093 float et;
0094 float hoe;
0095 } GCTtower_t;
0096
0097 typedef struct {
0098 GCTtower_t GCTtower[nBarrelEta / 2][nBarrelPhi];
0099 } GCTintTowers_t;
0100
0101 inline int getPeakBinOf3(float et0, float et1, float et2) {
0102 int x;
0103 float temp;
0104 if (et0 > et1) {
0105 x = 0;
0106 temp = et0;
0107 } else {
0108 x = 1;
0109 temp = et1;
0110 }
0111 if (et2 > temp) {
0112 x = 2;
0113 }
0114 return x;
0115 }
0116
0117 inline int getEtCenterOf3(float et0, float et1, float et2) {
0118 float etSum = et0 + et1 + et2;
0119 float iEtSum = 0.5 * et0 + 1.5 * et1 + 2.5 * et2;
0120 int iAve = 0xEEF;
0121 if (iEtSum <= etSum)
0122 iAve = 0;
0123 else if (iEtSum <= 2 * etSum)
0124 iAve = 1;
0125 else
0126 iAve = 2;
0127 return iAve;
0128 }
0129
0130 inline void makeST(const float GCTintTowers[nBarrelEta / 2][nBarrelPhi],
0131 GCTsupertower_t supertower_return[nSTEta][nSTPhi]) {
0132 float et_sumEta[nSTEta][nSTPhi][3];
0133 float stripEta[nSTEta][nSTPhi][3];
0134 float stripPhi[nSTEta][nSTPhi][3];
0135
0136 float ex_et[nBarrelEta / 2 + 1][nBarrelPhi];
0137 for (int j = 0; j < nBarrelPhi; j++) {
0138 ex_et[nBarrelEta / 2][j] = 0;
0139 for (int i = 0; i < nBarrelEta / 2; i++) {
0140 ex_et[i][j] = GCTintTowers[i][j];
0141 }
0142 }
0143
0144 int index_i = 0;
0145 int index_j = 0;
0146 for (int i = 0; i < nBarrelEta / 2 + 1; i += 3) {
0147 index_j = 0;
0148 for (int j = 0; j < nBarrelPhi; j += 3) {
0149 stripEta[index_i][index_j][0] = ex_et[i][j] + ex_et[i][j + 1] + ex_et[i][j + 2];
0150 stripEta[index_i][index_j][1] = ex_et[i + 1][j] + ex_et[i + 1][j + 1] + ex_et[i + 1][j + 2];
0151 stripEta[index_i][index_j][2] = ex_et[i + 2][j] + ex_et[i + 2][j + 1] + ex_et[i + 2][j + 2];
0152 stripPhi[index_i][index_j][0] = ex_et[i][j] + ex_et[i + 1][j] + ex_et[i + 2][j];
0153 stripPhi[index_i][index_j][1] = ex_et[i][j + 1] + ex_et[i + 1][j + 1] + ex_et[i + 2][j + 1];
0154 stripPhi[index_i][index_j][2] = ex_et[i][j + 2] + ex_et[i + 1][j + 2] + ex_et[i + 2][j + 2];
0155 for (int k = 0; k < 3; k++) {
0156 et_sumEta[index_i][index_j][k] = ex_et[i + k][j] + ex_et[i + k][j + 1] + ex_et[i + k][j + 2];
0157 }
0158 index_j++;
0159 }
0160 index_i++;
0161 }
0162 for (int i = 0; i < nSTEta; i++) {
0163 for (int j = 0; j < nSTPhi; j++) {
0164 GCTsupertower_t temp;
0165 temp.et = 0;
0166 temp.eta = 0;
0167 temp.phi = 0;
0168 temp.towerEta = 0;
0169 temp.towerPhi = 0;
0170 temp.towerEt = 0;
0171 temp.centerEta = 0;
0172 temp.centerPhi = 0;
0173 if (i < 6) {
0174 float supertowerEt = et_sumEta[i][j][0] + et_sumEta[i][j][1] + et_sumEta[i][j][2];
0175 temp.et = supertowerEt;
0176 temp.eta = i;
0177 temp.phi = j;
0178 int peakEta = getPeakBinOf3(stripEta[i][j][0], stripEta[i][j][1], stripEta[i][j][2]);
0179 temp.towerEta = peakEta;
0180 int peakPhi = getPeakBinOf3(stripPhi[i][j][0], stripPhi[i][j][1], stripPhi[i][j][2]);
0181 temp.towerPhi = peakPhi;
0182 float peakEt = ex_et[i * 3 + peakEta][j * 3 + peakPhi];
0183 temp.towerEt = peakEt;
0184 int cEta = getEtCenterOf3(stripEta[i][j][0], stripEta[i][j][1], stripEta[i][j][2]);
0185 temp.centerEta = cEta;
0186 int cPhi = getEtCenterOf3(stripPhi[i][j][0], stripPhi[i][j][1], stripPhi[i][j][2]);
0187 temp.centerPhi = cPhi;
0188 }
0189 supertower_return[i][j] = temp;
0190 }
0191 }
0192 }
0193
0194 inline void makeST_hgcal(const float hgcalTowers[nHgcalEta / 2][nHgcalPhi],
0195 GCTsupertower_t supertower_return[nSTEta][nSTPhi]) {
0196 float et_sumEta[nSTEta][nSTPhi][3];
0197 float stripEta[nSTEta][nSTPhi][3];
0198 float stripPhi[nSTEta][nSTPhi][3];
0199
0200 int index_i = 0;
0201 int index_j = 0;
0202 for (int i = 0; i < nHgcalEta / 2; i += 3) {
0203 index_j = 0;
0204 for (int j = 0; j < nHgcalPhi; j += 3) {
0205 stripEta[index_i][index_j][0] = hgcalTowers[i][j] + hgcalTowers[i][j + 1] + hgcalTowers[i][j + 2];
0206 stripEta[index_i][index_j][1] = hgcalTowers[i + 1][j] + hgcalTowers[i + 1][j + 1] + hgcalTowers[i + 1][j + 2];
0207 stripEta[index_i][index_j][2] = hgcalTowers[i + 2][j] + hgcalTowers[i + 2][j + 1] + hgcalTowers[i + 2][j + 2];
0208 stripPhi[index_i][index_j][0] = hgcalTowers[i][j] + hgcalTowers[i + 1][j] + hgcalTowers[i + 2][j];
0209 stripPhi[index_i][index_j][1] = hgcalTowers[i][j + 1] + hgcalTowers[i + 1][j + 1] + hgcalTowers[i + 2][j + 1];
0210 stripPhi[index_i][index_j][2] = hgcalTowers[i][j + 2] + hgcalTowers[i + 1][j + 2] + hgcalTowers[i + 2][j + 2];
0211 for (int k = 0; k < 3; k++) {
0212 et_sumEta[index_i][index_j][k] =
0213 hgcalTowers[i + k][j] + hgcalTowers[i + k][j + 1] + hgcalTowers[i + k][j + 2];
0214 }
0215 index_j++;
0216 }
0217 index_i++;
0218 }
0219
0220 for (int i = 0; i < nSTEta; i++) {
0221 for (int j = 0; j < nSTPhi; j++) {
0222 GCTsupertower_t temp;
0223 temp.et = 0;
0224 temp.eta = 0;
0225 temp.phi = 0;
0226 temp.towerEta = 0;
0227 temp.towerPhi = 0;
0228 temp.towerEt = 0;
0229 temp.centerEta = 0;
0230 temp.centerPhi = 0;
0231 if (i < 6) {
0232 float supertowerEt = et_sumEta[i][j][0] + et_sumEta[i][j][1] + et_sumEta[i][j][2];
0233 temp.et = supertowerEt;
0234 temp.eta = i;
0235 temp.phi = j;
0236 int peakEta = getPeakBinOf3(stripEta[i][j][0], stripEta[i][j][1], stripEta[i][j][2]);
0237 temp.towerEta = peakEta;
0238 int peakPhi = getPeakBinOf3(stripPhi[i][j][0], stripPhi[i][j][1], stripPhi[i][j][2]);
0239 temp.towerPhi = peakPhi;
0240 float peakEt = hgcalTowers[i * 3 + peakEta][j * 3 + peakPhi];
0241 temp.towerEt = peakEt;
0242 int cEta = getEtCenterOf3(stripEta[i][j][0], stripEta[i][j][1], stripEta[i][j][2]);
0243 temp.centerEta = cEta;
0244 int cPhi = getEtCenterOf3(stripPhi[i][j][0], stripPhi[i][j][1], stripPhi[i][j][2]);
0245 temp.centerPhi = cPhi;
0246 }
0247 supertower_return[i][j] = temp;
0248 }
0249 }
0250 }
0251
0252 inline void makeST_hf(const float hfTowers[nHfEta][nHfPhi], GCTsupertower_t supertower_return[nSTEta][nSTPhi]) {
0253 float et_sumEta[nSTEta][nSTPhi][3];
0254 float stripEta[nSTEta][nSTPhi][3];
0255 float stripPhi[nSTEta][nSTPhi][3];
0256
0257 int index_i = 0;
0258 int index_j = 0;
0259 for (int i = 0; i < nHfEta; i += 3) {
0260 index_j = 0;
0261 for (int j = 0; j < nHfPhi; j += 3) {
0262 stripEta[index_i][index_j][0] = hfTowers[i][j] + hfTowers[i][j + 1] + hfTowers[i][j + 2];
0263 stripEta[index_i][index_j][1] = hfTowers[i + 1][j] + hfTowers[i + 1][j + 1] + hfTowers[i + 1][j + 2];
0264 stripEta[index_i][index_j][2] = hfTowers[i + 2][j] + hfTowers[i + 2][j + 1] + hfTowers[i + 2][j + 2];
0265 stripPhi[index_i][index_j][0] = hfTowers[i][j] + hfTowers[i + 1][j] + hfTowers[i + 2][j];
0266 stripPhi[index_i][index_j][1] = hfTowers[i][j + 1] + hfTowers[i + 1][j + 1] + hfTowers[i + 2][j + 1];
0267 stripPhi[index_i][index_j][2] = hfTowers[i][j + 2] + hfTowers[i + 1][j + 2] + hfTowers[i + 2][j + 2];
0268 for (int k = 0; k < 3; k++) {
0269 et_sumEta[index_i][index_j][k] = hfTowers[i + k][j] + hfTowers[i + k][j + 1] + hfTowers[i + k][j + 2];
0270 }
0271 index_j++;
0272 }
0273 index_i++;
0274 }
0275
0276 for (int i = 0; i < nSTEta; i++) {
0277 for (int j = 0; j < nSTPhi; j++) {
0278 GCTsupertower_t temp;
0279 temp.et = 0;
0280 temp.eta = 0;
0281 temp.phi = 0;
0282 temp.towerEta = 0;
0283 temp.towerPhi = 0;
0284 temp.towerEt = 0;
0285 temp.centerEta = 0;
0286 temp.centerPhi = 0;
0287 float supertowerEt = et_sumEta[i][j][0] + et_sumEta[i][j][1] + et_sumEta[i][j][2];
0288 temp.et = supertowerEt;
0289 temp.eta = i;
0290 temp.phi = j;
0291 int peakEta = getPeakBinOf3(stripEta[i][j][0], stripEta[i][j][1], stripEta[i][j][2]);
0292 temp.towerEta = peakEta;
0293 int peakPhi = getPeakBinOf3(stripPhi[i][j][0], stripPhi[i][j][1], stripPhi[i][j][2]);
0294 temp.towerPhi = peakPhi;
0295 float peakEt = hfTowers[i * 3 + peakEta][j * 3 + peakPhi];
0296 temp.towerEt = peakEt;
0297 int cEta = getEtCenterOf3(stripEta[i][j][0], stripEta[i][j][1], stripEta[i][j][2]);
0298 temp.centerEta = cEta;
0299 int cPhi = getEtCenterOf3(stripPhi[i][j][0], stripPhi[i][j][1], stripPhi[i][j][2]);
0300 temp.centerPhi = cPhi;
0301 supertower_return[i][j] = temp;
0302 }
0303 }
0304 }
0305
0306 inline GCTsupertower_t bestOf2(const GCTsupertower_t& calotp0, const GCTsupertower_t& calotp1) {
0307 GCTsupertower_t x;
0308 x = (calotp0.et > calotp1.et) ? calotp0 : calotp1;
0309 return x;
0310 }
0311
0312 inline GCTsupertower_t getPeakBin24N(const etaStrip_t& etaStrip) {
0313 GCTsupertower_t best01 = bestOf2(etaStrip.cr[0], etaStrip.cr[1]);
0314 GCTsupertower_t best23 = bestOf2(etaStrip.cr[2], etaStrip.cr[3]);
0315 GCTsupertower_t best45 = bestOf2(etaStrip.cr[4], etaStrip.cr[5]);
0316 GCTsupertower_t best67 = bestOf2(etaStrip.cr[6], etaStrip.cr[7]);
0317 GCTsupertower_t best89 = bestOf2(etaStrip.cr[8], etaStrip.cr[9]);
0318 GCTsupertower_t best1011 = bestOf2(etaStrip.cr[10], etaStrip.cr[11]);
0319 GCTsupertower_t best1213 = bestOf2(etaStrip.cr[12], etaStrip.cr[13]);
0320 GCTsupertower_t best1415 = bestOf2(etaStrip.cr[14], etaStrip.cr[15]);
0321 GCTsupertower_t best1617 = bestOf2(etaStrip.cr[16], etaStrip.cr[17]);
0322 GCTsupertower_t best1819 = bestOf2(etaStrip.cr[18], etaStrip.cr[19]);
0323 GCTsupertower_t best2021 = bestOf2(etaStrip.cr[20], etaStrip.cr[21]);
0324 GCTsupertower_t best2223 = bestOf2(etaStrip.cr[22], etaStrip.cr[23]);
0325
0326 GCTsupertower_t best0123 = bestOf2(best01, best23);
0327 GCTsupertower_t best4567 = bestOf2(best45, best67);
0328 GCTsupertower_t best891011 = bestOf2(best89, best1011);
0329 GCTsupertower_t best12131415 = bestOf2(best1213, best1415);
0330 GCTsupertower_t best16171819 = bestOf2(best1617, best1819);
0331 GCTsupertower_t best20212223 = bestOf2(best2021, best2223);
0332
0333 GCTsupertower_t best01234567 = bestOf2(best0123, best4567);
0334 GCTsupertower_t best89101112131415 = bestOf2(best891011, best12131415);
0335 GCTsupertower_t best16to23 = bestOf2(best16171819, best20212223);
0336
0337 GCTsupertower_t best0to15 = bestOf2(best01234567, best89101112131415);
0338 GCTsupertower_t bestOf24 = bestOf2(best0to15, best16to23);
0339
0340 return bestOf24;
0341 }
0342
0343 inline towerMax getPeakBin8N(const etaStripPeak_t& etaStrip) {
0344 towerMax x;
0345
0346 GCTsupertower_t best01 = bestOf2(etaStrip.pk[0], etaStrip.pk[1]);
0347 GCTsupertower_t best23 = bestOf2(etaStrip.pk[2], etaStrip.pk[3]);
0348 GCTsupertower_t best45 = bestOf2(etaStrip.pk[4], etaStrip.pk[5]);
0349 GCTsupertower_t best67 = bestOf2(etaStrip.pk[6], etaStrip.pk[7]);
0350
0351 GCTsupertower_t best0123 = bestOf2(best01, best23);
0352 GCTsupertower_t best4567 = bestOf2(best45, best67);
0353 GCTsupertower_t bestOf8 = bestOf2(best0123, best4567);
0354
0355 x.energy = bestOf8.et;
0356 x.phi = bestOf8.phi;
0357 x.eta = bestOf8.eta;
0358 x.energyMax = bestOf8.towerEt;
0359 x.etaMax = bestOf8.towerEta;
0360 x.phiMax = bestOf8.towerPhi;
0361 x.etaCenter = bestOf8.centerEta;
0362 x.phiCenter = bestOf8.centerPhi;
0363 return x;
0364 }
0365
0366 inline jetInfo getJetPosition(GCTsupertower_t temp[nSTEta][nSTPhi]) {
0367 etaStripPeak_t etaStripPeak;
0368 jetInfo jet;
0369
0370 for (int i = 0; i < nSTEta; i++) {
0371 etaStrip_t test;
0372 for (int j = 0; j < nSTPhi; j++) {
0373 test.cr[j] = temp[i][j];
0374 }
0375 etaStripPeak.pk[i] = getPeakBin24N(test);
0376 }
0377
0378 towerMax peakIn8;
0379 peakIn8 = getPeakBin8N(etaStripPeak);
0380
0381 jet.seedEnergy = peakIn8.energy;
0382 jet.energy = 0;
0383 jet.tauEt = 0;
0384 jet.eta = peakIn8.eta;
0385 jet.phi = peakIn8.phi;
0386 jet.energyMax = peakIn8.energyMax;
0387 jet.etaMax = peakIn8.etaMax;
0388 jet.phiMax = peakIn8.phiMax;
0389 jet.etaCenter = peakIn8.etaCenter;
0390 jet.phiCenter = peakIn8.phiCenter;
0391
0392 return jet;
0393 }
0394
0395 inline jetInfo getJetValues(GCTsupertower_t tempX[nSTEta][nSTPhi], int seed_eta, int seed_phi) {
0396 float temp[nSTEta + 2][nSTPhi + 2];
0397 float eta_slice[3] = {0.};
0398 jetInfo jet_tmp;
0399
0400 for (int i = 0; i < nSTEta + 2; i++) {
0401 for (int k = 0; k < nSTPhi + 2; k++) {
0402 temp[i][k] = 0;
0403 }
0404 }
0405
0406 for (int i = 0; i < nSTEta; i++) {
0407 for (int k = 0; k < nSTPhi; k++) {
0408 temp[i + 1][k + 1] = tempX[i][k].et;
0409 }
0410 }
0411
0412 int seed_eta1, seed_phi1;
0413
0414 seed_eta1 = seed_eta;
0415 seed_phi1 = seed_phi;
0416 float tmp1, tmp2, tmp3;
0417
0418 for (int j = 0; j < nSTEta; j++) {
0419 for (int k = 0; k < nSTPhi; k++) {
0420 if (j == seed_eta1 && k == seed_phi1) {
0421 for (int m = 0; m < 3; m++) {
0422 tmp1 = temp[j + m][k];
0423 tmp2 = temp[j + m][k + 1];
0424 tmp3 = temp[j + m][k + 2];
0425 eta_slice[m] = tmp1 + tmp2 + tmp3;
0426 }
0427 }
0428 }
0429 }
0430
0431 jet_tmp.energy = eta_slice[0] + eta_slice[1] + eta_slice[2];
0432 jet_tmp.tauEt = eta_slice[1];
0433
0434 jet_tmp.etaCenter =
0435 3 * seed_eta + tempX[seed_eta][seed_phi].centerEta;
0436 jet_tmp.phiCenter =
0437 3 * seed_phi + tempX[seed_eta][seed_phi].centerPhi;
0438 jet_tmp.etaMax = 3 * seed_eta + tempX[seed_eta][seed_phi].towerEta;
0439 jet_tmp.phiMax = 3 * seed_phi + tempX[seed_eta][seed_phi].towerPhi;
0440
0441
0442 for (int i = 0; i < nSTEta; i++) {
0443 if (i + 1 >= seed_eta && i <= seed_eta + 1) {
0444 for (int k = 0; k < nSTPhi; k++) {
0445 if (k + 1 >= seed_phi && k <= seed_phi + 1)
0446 tempX[i][k].et = 0;
0447 }
0448 }
0449 }
0450
0451 return jet_tmp;
0452 }
0453
0454 inline jetInfo getRegion(GCTsupertower_t temp[nSTEta][nSTPhi], float TTseedThreshold) {
0455 jetInfo jet_tmp, jet;
0456 jet_tmp = getJetPosition(temp);
0457 int seed_phi = jet_tmp.phi;
0458 int seed_eta = jet_tmp.eta;
0459 float seed_tower_energy = jet_tmp.energyMax;
0460 jet = getJetValues(temp, seed_eta, seed_phi);
0461 if (seed_tower_energy > TTseedThreshold) {
0462 jet_tmp.energy = jet.energy;
0463 jet_tmp.tauEt = jet.tauEt;
0464 } else {
0465 jet_tmp.energy = 0.;
0466 jet_tmp.tauEt = 0.;
0467 }
0468 jet_tmp.etaCenter = jet.etaCenter;
0469 jet_tmp.phiCenter = jet.phiCenter;
0470 jet_tmp.etaMax = jet.etaMax;
0471 jet_tmp.phiMax = jet.phiMax;
0472 return jet_tmp;
0473 }
0474
0475 inline bool compareByEt(l1tp2::Phase2L1CaloJet i, l1tp2::Phase2L1CaloJet j) { return (i.jetEt() > j.jetEt()); };
0476
0477 }
0478
0479 #endif