File indexing completed on 2025-02-05 23:51:44
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "RecoMET/METProducers/interface/CaloMETProducer.h"
0010
0011 #include "FWCore/Framework/interface/ConsumesCollector.h"
0012
0013 #include "FWCore/Utilities/interface/InputTag.h"
0014
0015 #include "DataFormats/Common/interface/Handle.h"
0016
0017 #include "DataFormats/METReco/interface/METFwd.h"
0018 #include "DataFormats/METReco/interface/CaloMETFwd.h"
0019 #include "DataFormats/METReco/interface/CaloMET.h"
0020 #include "DataFormats/METReco/interface/CommonMETData.h"
0021
0022 #include "RecoMET/METAlgorithms/interface/METAlgo.h"
0023 #include "RecoMET/METAlgorithms/interface/SignAlgoResolutions.h"
0024 #include "RecoMET/METAlgorithms/interface/CaloSpecificAlgo.h"
0025 #include "RecoMET/METAlgorithms/interface/SignCaloSpecificAlgo.h"
0026
0027 #include <cstring>
0028
0029
0030 namespace cms {
0031
0032
0033 CaloMETProducer::CaloMETProducer(const edm::ParameterSet& iConfig)
0034 : inputToken_(consumes<edm::View<reco::Candidate> >(iConfig.getParameter<edm::InputTag>("src"))),
0035 calculateSignificance_(iConfig.getParameter<bool>("calculateSignificance")),
0036 resolutions_(nullptr),
0037 globalThreshold_(iConfig.getParameter<double>("globalThreshold")) {
0038 noHF_ = iConfig.getParameter<bool>("noHF");
0039
0040 std::string alias = iConfig.getParameter<std::string>("alias");
0041 produces<reco::CaloMETCollection>().setBranchAlias(alias);
0042
0043 if (calculateSignificance_)
0044 resolutions_ = new metsig::SignAlgoResolutions(iConfig);
0045 }
0046
0047
0048 CaloMETProducer::~CaloMETProducer() {
0049 if (resolutions_)
0050 delete resolutions_;
0051 }
0052
0053
0054 void CaloMETProducer::produce(edm::Event& event, const edm::EventSetup& setup) {
0055 edm::Handle<edm::View<reco::Candidate> > input;
0056 event.getByToken(inputToken_, input);
0057
0058 METAlgo algo;
0059 CommonMETData commonMETdata = algo.run(*input.product(), globalThreshold_);
0060
0061 CaloSpecificAlgo calospecalgo;
0062 reco::CaloMET calomet = calospecalgo.addInfo(input, commonMETdata, noHF_, globalThreshold_);
0063
0064
0065
0066
0067
0068
0069
0070
0071
0072
0073 auto calometcoll = std::make_unique<reco::CaloMETCollection>();
0074 calometcoll->push_back(calomet);
0075 event.put(std::move(calometcoll));
0076 }
0077
0078
0079 void CaloMETProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0080 edm::ParameterSetDescription desc;
0081 desc.add<edm::InputTag>("src", edm::InputTag("towerMaker"));
0082 desc.add<bool>("calculateSignificance", false);
0083 desc.add<double>("globalThreshold", 0.3);
0084 desc.add<bool>("noHF", false);
0085 desc.add<std::string>("alias", "");
0086 metsig::SignAlgoResolutions::fillPSetDescription(desc);
0087 descriptions.addWithDefaultLabel(desc);
0088 }
0089
0090
0091 DEFINE_FWK_MODULE(CaloMETProducer);
0092 }
0093
0094