File indexing completed on 2023-03-17 11:12:19
0001
0002 #include "FWCore/Framework/interface/Frameworkfwd.h"
0003 #include "L1Trigger/L1THGCal/interface/backend/HGCalTriggerClusterInterpreterBase.h"
0004 #include "DataFormats/L1THGCal/interface/HGCalMulticluster.h"
0005 #include "L1Trigger/L1THGCal/interface/HGCalTriggerTools.h"
0006
0007 class HGCalTriggerClusterInterpretationEM : public HGCalTriggerClusterInterpreterBase {
0008 public:
0009 HGCalTriggerClusterInterpretationEM();
0010 ~HGCalTriggerClusterInterpretationEM() override{};
0011 void initialize(const edm::ParameterSet& conf) final;
0012 void setGeometry(const HGCalTriggerGeometryBase* const) final;
0013 void interpret(l1t::HGCalMulticlusterBxCollection& multiclusters) const final;
0014
0015 private:
0016 std::vector<double> layer_containment_corrs_;
0017 std::vector<double> scale_corrections_coeff_;
0018 std::vector<double> dr_bylayer_;
0019
0020 HGCalTriggerTools triggerTools_;
0021 };
0022
0023 DEFINE_HGC_TPG_CLUSTER_INTERPRETER(HGCalTriggerClusterInterpretationEM, "HGCalTriggerClusterInterpretationEM");
0024
0025 HGCalTriggerClusterInterpretationEM::HGCalTriggerClusterInterpretationEM() {}
0026
0027 void HGCalTriggerClusterInterpretationEM::initialize(const edm::ParameterSet& conf) {
0028 layer_containment_corrs_ = conf.getParameter<std::vector<double>>("layer_containment_corrs");
0029 scale_corrections_coeff_ = conf.getParameter<std::vector<double>>("scale_correction_coeff");
0030 dr_bylayer_ = conf.getParameter<std::vector<double>>("dr_bylayer");
0031
0032 const unsigned corrections_size = 2;
0033 if (scale_corrections_coeff_.size() != corrections_size) {
0034 throw cms::Exception("HGCTriggerParameterError")
0035 << "HGCalTriggerClusterInterpretationEM::scale_correction_coeff parameter has size: "
0036 << scale_corrections_coeff_.size() << " while expected is " << corrections_size;
0037 }
0038 if (layer_containment_corrs_.size() != dr_bylayer_.size()) {
0039 throw cms::Exception("HGCTriggerParameterError")
0040 << "HGCalTriggerClusterInterpretationEM::layer_containment_corrs and "
0041 "HGCalTriggerClusterInterpretationEM::dr_bylayer have different size!";
0042 }
0043 }
0044
0045 void HGCalTriggerClusterInterpretationEM::setGeometry(const HGCalTriggerGeometryBase* const geom) {
0046 triggerTools_.setGeometry(geom);
0047 }
0048
0049 void HGCalTriggerClusterInterpretationEM::interpret(l1t::HGCalMulticlusterBxCollection& multiclusters) const {
0050 for (unsigned int idx = 0; idx != multiclusters.size(); idx++) {
0051 l1t::HGCalMulticluster& cluster3d = multiclusters[idx];
0052
0053 const GlobalPoint& cluster3d_position = cluster3d.centreProj();
0054 double energy = 0.;
0055
0056 for (const auto& cluster2d : cluster3d.constituents()) {
0057 const unsigned layer = triggerTools_.triggerLayer(cluster2d.first);
0058 if (layer <= layer_containment_corrs_.size() - 1) {
0059 double dr = (cluster3d_position - cluster2d.second->centreProj()).mag();
0060 if (dr <= dr_bylayer_.at(layer)) {
0061 energy += layer_containment_corrs_.at(layer) * cluster2d.second->energy();
0062 }
0063 }
0064 }
0065 energy += scale_corrections_coeff_.at(1) * std::abs(cluster3d.eta()) + scale_corrections_coeff_.at(0);
0066 cluster3d.saveEnergyInterpretation(l1t::HGCalMulticluster::EnergyInterpretation::EM, max(energy, 0.));
0067 }
0068 }