File indexing completed on 2024-07-03 04:18:02
0001 #ifndef L1Trigger_Phase2L1GT_L1GTScales_h
0002 #define L1Trigger_Phase2L1GT_L1GTScales_h
0003
0004 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006
0007 #include <cmath>
0008
0009 namespace l1t {
0010 class L1GTScales {
0011 public:
0012 static constexpr int RELATIVE_ISOLATION_RESOLUTION = 18;
0013
0014
0015 static constexpr int INV_MASS_SQR_OVER_2_DR_SQR_RESOLUTION = 19;
0016
0017 L1GTScales(double pT_lsb,
0018 double phi_lsb,
0019 double eta_lsb,
0020 double z0_lsb,
0021
0022 double isolationPT_lsb,
0023 double beta_lsb,
0024 double mass_lsb,
0025 double seed_pT_lsb,
0026 double seed_dZ_lsb,
0027 double scalarSumPT_lsb,
0028 double sum_pT_pv_lsb,
0029 int pos_chg,
0030 int neg_chg);
0031
0032 L1GTScales(const edm::ParameterSet &);
0033
0034 static void fillPSetDescription(edm::ParameterSetDescription &);
0035
0036 int to_hw_pT_ceil(double value) const { return std::ceil(value / pT_lsb_); };
0037 int to_hw_phi_ceil(double value) const { return std::ceil(value / phi_lsb_); };
0038 int to_hw_eta_ceil(double value) const { return std::ceil(value / eta_lsb_); };
0039 int to_hw_z0_ceil(double value) const { return std::ceil(value / z0_lsb_); };
0040
0041 int to_hw_relative_isolationPT_ceil(double value) const {
0042 return std::ceil(pT_lsb_ * value * std::pow(2, RELATIVE_ISOLATION_RESOLUTION) / isolationPT_lsb_);
0043 }
0044 int to_hw_isolationPT_ceil(double value) const { return std::ceil(value / isolationPT_lsb_); }
0045 int to_hw_beta_ceil(double value) const { return std::ceil(value / beta_lsb_); };
0046 int to_hw_mass_ceil(double value) const { return std::ceil(value / mass_lsb_); };
0047 int to_hw_seed_pT_ceil(double value) const { return std::ceil(value / seed_pT_lsb_); };
0048 int to_hw_seed_z0_ceil(double value) const { return std::ceil(value / seed_z0_lsb_); };
0049 int to_hw_scalarSumPT_ceil(double value) const { return std::ceil(value / scalarSumPT_lsb_); };
0050 int to_hw_sum_pT_pv_ceil(double value) const { return std::ceil(value / sum_pT_pv_lsb_); };
0051
0052 int to_hw_dRSquared_ceil(double value) const { return std::ceil(value * value / (eta_lsb_ * eta_lsb_)); }
0053
0054 double to_hw_InvMassSqrDiv2(double value) const { return value * value / (2 * pT_lsb_ * pT_lsb_); }
0055 double to_hw_TransMassSqrDiv2(double value) const { return value * value / (2 * pT_lsb_ * pT_lsb_); }
0056
0057 double to_hw_PtSquared(double value) const { return value * value / (pT_lsb_ * pT_lsb_); }
0058 double to_hw_InvMassSqrOver2DR(double value) const {
0059 return value * value * eta_lsb_ * eta_lsb_ * std::pow(2, INV_MASS_SQR_OVER_2_DR_SQR_RESOLUTION) /
0060 (2 * pT_lsb_ * pT_lsb_);
0061 }
0062
0063 int to_hw_pT_floor(double value) const { return std::floor(value / pT_lsb_); };
0064 int to_hw_phi_floor(double value) const { return std::floor(value / phi_lsb_); };
0065 int to_hw_eta_floor(double value) const { return std::floor(value / eta_lsb_); };
0066 int to_hw_z0_floor(double value) const { return std::floor(value / z0_lsb_); };
0067
0068 int to_hw_relative_isolationPT_floor(double value) const {
0069 return std::floor(pT_lsb_ * value * std::pow(2, RELATIVE_ISOLATION_RESOLUTION) / isolationPT_lsb_);
0070 }
0071 int to_hw_isolationPT_floor(double value) const { return std::floor(value / isolationPT_lsb_); }
0072 int to_hw_beta_floor(double value) const { return std::floor(value / beta_lsb_); };
0073 int to_hw_mass_floor(double value) const { return std::floor(value / mass_lsb_); };
0074 int to_hw_seed_pT_floor(double value) const { return std::floor(value / seed_pT_lsb_); };
0075 int to_hw_seed_z0_floor(double value) const { return std::floor(value / seed_z0_lsb_); };
0076 int to_hw_scalarSumPT_floor(double value) const { return std::floor(value / scalarSumPT_lsb_); };
0077 int to_hw_sum_pT_pv_floor(double value) const { return std::floor(value / sum_pT_pv_lsb_); };
0078
0079 int to_hw_dRSquared_floor(double value) const { return std::floor(value * value / (eta_lsb_ * eta_lsb_)); }
0080
0081 double to_pT(int value) const { return value * pT_lsb_; };
0082 double to_phi(int value) const { return value * phi_lsb_; };
0083 double to_eta(int value) const { return value * eta_lsb_; };
0084 double to_z0(int value) const { return value * z0_lsb_; };
0085 double to_isolationPT(int value) const { return value * isolationPT_lsb_; }
0086 double to_scalarSumPT(int value) const { return value * scalarSumPT_lsb_; };
0087 int to_chg(int value) const { return value == pos_chg_ ? +1 : value == neg_chg_ ? -1 : 0; }
0088
0089 double pT_lsb() const { return pT_lsb_; }
0090 double phi_lsb() const { return phi_lsb_; }
0091 double eta_lsb() const { return eta_lsb_; }
0092 double z0_lsb() const { return z0_lsb_; }
0093 double isolationPT_lsb() const { return isolationPT_lsb_; }
0094
0095 double beta_lsb() const { return beta_lsb_; }
0096 double mass_lsb() const { return mass_lsb_; }
0097 double seed_pT_lsb() const { return seed_pT_lsb_; }
0098 double seed_z0_lsb() const { return seed_z0_lsb_; }
0099 double scalarSumPT_lsb() const { return scalarSumPT_lsb_; }
0100 double sum_pT_pv_lsb() const { return sum_pT_pv_lsb_; }
0101 int pos_chg() const { return pos_chg_; }
0102 int neg_chg() const { return neg_chg_; }
0103
0104 private:
0105 const double pT_lsb_;
0106 const double phi_lsb_;
0107 const double eta_lsb_;
0108 const double z0_lsb_;
0109
0110 const double isolationPT_lsb_;
0111 const double beta_lsb_;
0112 const double mass_lsb_;
0113 const double seed_pT_lsb_;
0114 const double seed_z0_lsb_;
0115 const double scalarSumPT_lsb_;
0116 const double sum_pT_pv_lsb_;
0117 const int pos_chg_;
0118 const int neg_chg_;
0119 };
0120 }
0121
0122 #endif