File indexing completed on 2024-04-06 12:24:11
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 #include "FWCore/Framework/interface/global/EDProducer.h"
0013 #include "FWCore/Framework/interface/Event.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/Utilities/interface/InputTag.h"
0016
0017 #include "DataFormats/Common/interface/ValueMap.h"
0018 #include "DataFormats/Common/interface/View.h"
0019
0020 #include "DataFormats/Candidate/interface/CandidateFwd.h"
0021 #include "DataFormats/Candidate/interface/Candidate.h"
0022 #include "PhysicsTools/TagAndProbe/interface/ColinsSoperVariables.h"
0023 #include "TLorentzVector.h"
0024
0025 class ColinsSoperVariablesComputer : public edm::global::EDProducer<> {
0026 public:
0027 explicit ColinsSoperVariablesComputer(const edm::ParameterSet& iConfig);
0028 ~ColinsSoperVariablesComputer() override;
0029
0030 void produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const override;
0031
0032 private:
0033 const edm::EDGetTokenT<edm::View<reco::Candidate>> parentBosonToken_;
0034 };
0035
0036 ColinsSoperVariablesComputer::ColinsSoperVariablesComputer(const edm::ParameterSet& iConfig)
0037 : parentBosonToken_(consumes<edm::View<reco::Candidate>>(iConfig.getParameter<edm::InputTag>("parentBoson"))) {
0038 produces<edm::ValueMap<float>>("costheta");
0039 produces<edm::ValueMap<float>>("sin2theta");
0040 produces<edm::ValueMap<float>>("tanphi");
0041 }
0042
0043 ColinsSoperVariablesComputer::~ColinsSoperVariablesComputer() {}
0044
0045 void ColinsSoperVariablesComputer::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0046 using namespace edm;
0047
0048
0049 Handle<View<reco::Candidate>> bosons;
0050 iEvent.getByToken(parentBosonToken_, bosons);
0051
0052
0053 std::vector<float> values;
0054 std::vector<float> values2;
0055 std::vector<float> values3;
0056
0057
0058 double costheta = -10.0;
0059 double sin2theta = -10.0;
0060 double tanphi = -10.0;
0061
0062 const reco::Candidate* daughter1 = nullptr;
0063 const reco::Candidate* daughter2 = nullptr;
0064 TLorentzVector mu(0., 0., 0., 0.);
0065 TLorentzVector mubar(0., 0., 0., 0.);
0066 bool isOS = false;
0067 int charge1 = 0, charge2 = 0;
0068 double res[3] = {-10., -10., -10.};
0069
0070 View<reco::Candidate>::const_iterator boson, endbosons = bosons->end();
0071
0072 for (boson = bosons->begin(); boson != endbosons; ++boson) {
0073 daughter1 = boson->daughter(0);
0074 daughter2 = boson->daughter(1);
0075
0076 if (!(nullptr == daughter1 || nullptr == daughter2)) {
0077 charge1 = daughter1->charge();
0078 charge2 = daughter2->charge();
0079 isOS = charge1 * charge2 < 0;
0080 if (isOS && charge1 < 0) {
0081 mu.SetPxPyPzE(daughter1->px(), daughter1->py(), daughter1->pz(), daughter1->energy());
0082 mubar.SetPxPyPzE(daughter2->px(), daughter2->py(), daughter2->pz(), daughter2->energy());
0083 }
0084 if (isOS && charge1 > 0) {
0085 mu.SetPxPyPzE(daughter2->px(), daughter2->py(), daughter2->pz(), daughter2->energy());
0086 mubar.SetPxPyPzE(daughter1->px(), daughter1->py(), daughter1->pz(), daughter1->energy());
0087 }
0088 }
0089
0090 calCSVariables(mu, mubar, res, boson->pz() < 0.0);
0091
0092 costheta = res[0];
0093 sin2theta = res[1];
0094 tanphi = res[2];
0095
0096 values.push_back(costheta);
0097 values2.push_back(sin2theta);
0098 values3.push_back(tanphi);
0099 }
0100
0101
0102 auto valMap = std::make_unique<ValueMap<float>>();
0103 ValueMap<float>::Filler filler(*valMap);
0104 filler.insert(bosons, values.begin(), values.end());
0105 filler.fill();
0106 iEvent.put(std::move(valMap), "costheta");
0107
0108
0109 auto valMap2 = std::make_unique<ValueMap<float>>();
0110 ValueMap<float>::Filler filler2(*valMap2);
0111 filler2.insert(bosons, values2.begin(), values2.end());
0112 filler2.fill();
0113 iEvent.put(std::move(valMap2), "sin2theta");
0114
0115
0116 auto valMap3 = std::make_unique<ValueMap<float>>();
0117 ValueMap<float>::Filler filler3(*valMap3);
0118 filler3.insert(bosons, values3.begin(), values3.end());
0119 filler3.fill();
0120 iEvent.put(std::move(valMap3), "tanphi");
0121 }
0122
0123 #include "FWCore/Framework/interface/MakerMacros.h"
0124 DEFINE_FWK_MODULE(ColinsSoperVariablesComputer);