File indexing completed on 2021-03-25 23:59:57
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "FWCore/Framework/interface/Event.h"
0015 #include "FWCore/Framework/interface/EventSetup.h"
0016 #include "FWCore/Framework/interface/global/EDFilter.h"
0017
0018 #include "DataFormats/CTPPSDetId/interface/CTPPSDetId.h"
0019 #include "DataFormats/CTPPSReco/interface/CTPPSLocalTrackLite.h"
0020 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0021 #include "DataFormats/ProtonReco/interface/ForwardProton.h"
0022
0023 #include "CondFormats/RunInfo/interface/LHCInfo.h"
0024 #include "CondFormats/DataRecord/interface/LHCInfoRcd.h"
0025
0026 #include "FWCore/Framework/interface/MakerMacros.h"
0027
0028
0029
0030 class HLTPPSJetComparisonFilter : public edm::global::EDFilter<> {
0031 public:
0032 explicit HLTPPSJetComparisonFilter(const edm::ParameterSet &);
0033 ~HLTPPSJetComparisonFilter() override;
0034
0035 static void fillDescriptions(edm::ConfigurationDescriptions &);
0036 bool filter(edm::StreamID, edm::Event &, const edm::EventSetup &) const override;
0037
0038 private:
0039
0040 edm::ESGetToken<LHCInfo, LHCInfoRcd> const lhcInfoRcdToken_;
0041
0042 edm::ParameterSet param_;
0043
0044 edm::InputTag jetInputTag_;
0045 edm::EDGetTokenT<reco::PFJetCollection> jet_token_;
0046
0047 edm::InputTag forwardProtonInputTag_;
0048 edm::EDGetTokenT<std::vector<reco::ForwardProton>> recoProtonSingleRPToken_;
0049
0050 std::string lhcInfoLabel_;
0051
0052 double maxDiffxi_;
0053 double maxDiffm_;
0054 double maxDiffy_;
0055
0056 unsigned int n_jets_;
0057
0058 bool do_xi_;
0059 bool do_my_;
0060 };
0061
0062
0063
0064 void HLTPPSJetComparisonFilter::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0065 edm::ParameterSetDescription desc;
0066
0067 desc.add<edm::InputTag>("jetInputTag", edm::InputTag("hltAK4PFJetsCorrected"))
0068 ->setComment("input tag of the jet track collection");
0069 desc.add<edm::InputTag>("forwardProtonInputTag", edm::InputTag("ctppsProtons", "singleRP"))
0070 ->setComment("input tag of the forward proton collection");
0071
0072 desc.add<std::string>("lhcInfoLabel", std::string(""))->setComment("label used for LHCInfo");
0073
0074 desc.add<double>("maxDiffxi", 1.)
0075 ->setComment("maximum relative deviation of RP xi from dijet xi. Used with do_xi option");
0076 desc.add<double>("maxDiffm", 1.)
0077 ->setComment("maximum relative deviation of RP m from dijet m- Used with do_my option");
0078 desc.add<double>("maxDiffy", 1.)
0079 ->setComment("maximum absolute deviation of RP y from dijet y. Used with do_my option");
0080
0081 desc.add<unsigned int>("nJets", 2)->setComment("number of jets to be used");
0082
0083 desc.add<bool>("do_xi", true)->setComment("flag to require xi matching");
0084 desc.add<bool>("do_my", false)->setComment("flag to require m,y matching");
0085
0086 descriptions.addWithDefaultLabel(desc);
0087 return;
0088 }
0089
0090
0091
0092 HLTPPSJetComparisonFilter::~HLTPPSJetComparisonFilter() = default;
0093
0094 HLTPPSJetComparisonFilter::HLTPPSJetComparisonFilter(const edm::ParameterSet &iConfig)
0095 : lhcInfoRcdToken_(esConsumes()),
0096 jetInputTag_(iConfig.getParameter<edm::InputTag>("jetInputTag")),
0097 jet_token_(consumes<reco::PFJetCollection>(jetInputTag_)),
0098
0099 forwardProtonInputTag_(iConfig.getParameter<edm::InputTag>("forwardProtonInputTag")),
0100 recoProtonSingleRPToken_(consumes<std::vector<reco::ForwardProton>>(forwardProtonInputTag_)),
0101
0102 lhcInfoLabel_(iConfig.getParameter<std::string>("lhcInfoLabel")),
0103
0104 maxDiffxi_(iConfig.getParameter<double>("maxDiffxi")),
0105 maxDiffm_(iConfig.getParameter<double>("maxDiffm")),
0106 maxDiffy_(iConfig.getParameter<double>("maxDiffy")),
0107
0108 n_jets_(iConfig.getParameter<unsigned int>("nJets")),
0109
0110 do_xi_(iConfig.getParameter<bool>("do_xi")),
0111 do_my_(iConfig.getParameter<bool>("do_my")) {}
0112
0113
0114
0115 bool HLTPPSJetComparisonFilter::filter(edm::StreamID, edm::Event &iEvent, const edm::EventSetup &iSetup) const {
0116 auto const &hLHCInfo = iSetup.getHandle(lhcInfoRcdToken_);
0117 float sqs = 2. * hLHCInfo->energy();
0118
0119 edm::Handle<reco::PFJetCollection> jets;
0120 iEvent.getByToken(jet_token_, jets);
0121
0122 edm::Handle<std::vector<reco::ForwardProton>> recoSingleRPProtons;
0123 iEvent.getByToken(recoProtonSingleRPToken_, recoSingleRPProtons);
0124
0125 if (jets->size() < n_jets_)
0126 return false;
0127
0128 if (do_xi_ && maxDiffxi_ > 0) {
0129
0130 float sum45 = 0, sum56 = 0;
0131
0132 for (unsigned int i = 0; i < n_jets_; i++) {
0133 sum45 += (*jets)[i].energy() + (*jets)[i].pz();
0134 sum56 += (*jets)[i].energy() - (*jets)[i].pz();
0135 }
0136
0137 const float xi45 = sum45 / sqs;
0138 const float xi56 = sum56 / sqs;
0139
0140 float min45 = 1000., min56 = 1000.;
0141
0142 for (const auto &proton : *recoSingleRPProtons)
0143 {
0144 if (proton.validFit())
0145 {
0146 const auto &xi = proton.xi();
0147
0148 CTPPSDetId rpId(
0149 (*proton.contributingLocalTracks().begin())->rpId());
0150
0151 if (rpId.arm() == 0 && std::abs(xi - xi45) < min45)
0152 min45 = std::abs(xi - xi45);
0153 if (rpId.arm() == 1 && std::abs(xi - xi56) < min56)
0154 min56 = std::abs(xi - xi56);
0155 }
0156 }
0157
0158 if (min56 / xi56 > maxDiffxi_ || min45 / xi45 > maxDiffxi_)
0159 return false;
0160 }
0161
0162 if (do_my_) {
0163
0164
0165 ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<float>> j_sum;
0166 for (unsigned int i = 0; i < n_jets_; i++)
0167 j_sum = j_sum + (*jets)[i].p4();
0168
0169 const auto &mjet = j_sum.M();
0170 const auto &yjet = j_sum.Rapidity();
0171
0172 for (const auto &proton1 : *recoSingleRPProtons)
0173 {
0174 if (proton1.validFit()) {
0175 CTPPSDetId rpId1(
0176 (*proton1.contributingLocalTracks().begin())->rpId());
0177 if (rpId1.arm() == 0) {
0178 const auto &xi_45 = proton1.xi();
0179
0180 for (const auto &proton2 : *recoSingleRPProtons)
0181 {
0182 if (proton2.validFit()) {
0183 CTPPSDetId rpId2((*proton2.contributingLocalTracks().begin())->rpId());
0184 if (rpId2.arm() == 1) {
0185 const auto &xi_56 = proton2.xi();
0186
0187
0188 const auto &m = sqs * sqrt(xi_45 * xi_56);
0189 const auto &y = 0.5 * log(xi_45 / xi_56);
0190 if ((std::abs(m - mjet) / mjet < maxDiffm_ || maxDiffm_ <= 0) &&
0191 (std::abs(y - yjet) < maxDiffy_ || maxDiffy_ <= 0))
0192 return true;
0193 }
0194 }
0195 }
0196 }
0197 }
0198 }
0199 return false;
0200 }
0201
0202 return true;
0203 }
0204
0205 DEFINE_FWK_MODULE(HLTPPSJetComparisonFilter);