File indexing completed on 2022-02-18 08:23:38
0001 #include "L1Trigger/L1THGCal/interface/backend/HGCalStage1TruncationImpl_SA.h"
0002 #include <cmath>
0003
0004 unsigned HGCalStage1TruncationImplSA::run(const l1thgcfirmware::HGCalTriggerCellSACollection& tcs_in,
0005 const l1thgcfirmware::Stage1TruncationConfig& theConf,
0006 l1thgcfirmware::HGCalTriggerCellSACollection& tcs_out) const {
0007 unsigned sector120 = theConf.phiSector();
0008 std::unordered_map<unsigned, l1thgcfirmware::HGCalTriggerCellSACollection> tcs_per_bin;
0009
0010
0011 bool do_truncate = theConf.doTruncate();
0012 double rozmin = theConf.rozMin();
0013 double rozmax = theConf.rozMax();
0014 unsigned rozbins = theConf.rozBins();
0015 const std::vector<unsigned>& maxtcsperbin = theConf.maxTcsPerBin();
0016 const std::vector<double>& phiedges = theConf.phiEdges();
0017
0018 constexpr double margin = 1.001;
0019 double roz_bin_size = (rozbins > 0 ? (rozmax - rozmin) * margin / double(rozbins) : 0.);
0020
0021
0022 for (const auto& tc : tcs_in) {
0023 double x = tc.x();
0024 double y = tc.y();
0025 double z = tc.z();
0026 unsigned roverzbin = 0;
0027 if (roz_bin_size > 0.) {
0028 double roverz = std::sqrt(x * x + y * y) / std::abs(z) - rozmin;
0029 roverz = std::clamp(roverz, 0., rozmax - rozmin);
0030 roverzbin = unsigned(roverz / roz_bin_size);
0031 }
0032 double phi = rotatedphi(x, y, z, sector120);
0033 int phibin = phiBin(roverzbin, phi, phiedges);
0034 if (phibin < 0)
0035 return 1;
0036 unsigned packed_bin = packBin(roverzbin, phibin);
0037
0038 tcs_per_bin[packed_bin].push_back(tc);
0039 }
0040
0041 for (auto& bin_tcs : tcs_per_bin) {
0042 std::sort(bin_tcs.second.begin(),
0043 bin_tcs.second.end(),
0044 [](const l1thgcfirmware::HGCalTriggerCell& a, const l1thgcfirmware::HGCalTriggerCell& b) -> bool {
0045 return a.mipPt() > b.mipPt();
0046 });
0047
0048 unsigned roverzbin = 0;
0049 unsigned phibin = 0;
0050 unpackBin(bin_tcs.first, roverzbin, phibin);
0051 if (roverzbin >= maxtcsperbin.size())
0052 return 1;
0053
0054 unsigned max_tc = maxtcsperbin[roverzbin];
0055 if (do_truncate && bin_tcs.second.size() > max_tc) {
0056 bin_tcs.second.resize(max_tc);
0057 }
0058
0059 for (const auto& tc : bin_tcs.second) {
0060 tcs_out.push_back(tc);
0061 }
0062 }
0063
0064 return 0;
0065 }
0066
0067 unsigned HGCalStage1TruncationImplSA::packBin(unsigned roverzbin, unsigned phibin) const {
0068 unsigned packed_bin = 0;
0069 packed_bin |= ((roverzbin & mask_roz_) << offset_roz_);
0070 packed_bin |= (phibin & mask_phi_);
0071 return packed_bin;
0072 }
0073
0074 void HGCalStage1TruncationImplSA::unpackBin(unsigned packedbin, unsigned& roverzbin, unsigned& phibin) const {
0075 roverzbin = ((packedbin >> offset_roz_) & mask_roz_);
0076 phibin = (packedbin & mask_phi_);
0077 }
0078
0079 int HGCalStage1TruncationImplSA::phiBin(unsigned roverzbin, double phi, const std::vector<double>& phiedges) const {
0080 int phi_bin = 0;
0081 if (roverzbin >= phiedges.size())
0082 return -1;
0083 double phi_edge = phiedges[roverzbin];
0084 if (phi > phi_edge)
0085 phi_bin = 1;
0086 return phi_bin;
0087 }
0088
0089 double HGCalStage1TruncationImplSA::rotatedphi(double x, double y, double z, int sector) const {
0090 if (z > 0)
0091 x = -x;
0092 double phi = std::atan2(y, x);
0093
0094 if (sector == 1) {
0095 if (phi < M_PI and phi > 0)
0096 phi = phi - (2. * M_PI / 3.);
0097 else
0098 phi = phi + (4. * M_PI / 3.);
0099 } else if (sector == 2) {
0100 phi = phi + (2. * M_PI / 3.);
0101 }
0102 return phi;
0103 }