File indexing completed on 2023-10-25 09:55:41
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014 #include <memory>
0015
0016
0017 #include "FWCore/Framework/interface/Frameworkfwd.h"
0018 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0019 #include "FWCore/Framework/interface/Event.h"
0020 #include "FWCore/Framework/interface/EventSetup.h"
0021 #include "FWCore/Framework/interface/MakerMacros.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0024 #include "FWCore/Framework/interface/ESHandle.h"
0025
0026
0027 #include "JetMETCorrections/JetCorrector/interface/JetCorrector.h"
0028 #include "CondFormats/EcalObjects/interface/EcalChannelStatus.h"
0029 #include "CondFormats/DataRecord/interface/EcalChannelStatusRcd.h"
0030
0031
0032 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
0033 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0034 #include "DataFormats/JetReco/interface/JetID.h"
0035 #include "DataFormats/METReco/interface/PFMETCollection.h"
0036 #include "DataFormats/MuonReco/interface/Muon.h"
0037 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0038 #include "DataFormats/METReco/interface/PFMET.h"
0039 #include "DataFormats/MuonReco/interface/Muon.h"
0040 #include "DataFormats/METReco/interface/CaloMETCollection.h"
0041 #include "DataFormats/METReco/interface/CaloMET.h"
0042 #include "DataFormats/PatCandidates/interface/Jet.h"
0043
0044
0045 #include "FWCore/ServiceRegistry/interface/Service.h"
0046 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0047 #include "TH1.h"
0048 #include "TTree.h"
0049 #include "TF1.h"
0050 #include <TVector2.h>
0051
0052
0053 #include "L1Trigger/L1TNtuples/interface/L1AnalysisRecoJetDataFormat.h"
0054 #include "L1Trigger/L1TNtuples/interface/L1AnalysisRecoMetDataFormat.h"
0055
0056
0057
0058
0059
0060 class L1JetRecoTreeProducer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0061 public:
0062 explicit L1JetRecoTreeProducer(const edm::ParameterSet&);
0063 ~L1JetRecoTreeProducer() override;
0064
0065 private:
0066 void beginJob(void) override;
0067 void analyze(const edm::Event&, const edm::EventSetup&) override;
0068 void endJob() override;
0069
0070 void doPFJets(edm::Handle<reco::PFJetCollection> pfJets);
0071 void doPFJetCorr(edm::Handle<reco::PFJetCollection> pfJets, edm::Handle<reco::JetCorrector> pfJetCorr);
0072 void doPUPPIJets(edm::Handle<reco::PFJetCollection> puppiJets);
0073 void doCorrPUPPIJets(edm::Handle<std::vector<pat::Jet> > corrPuppiJets);
0074 void doCaloJets(edm::Handle<reco::CaloJetCollection> caloJets);
0075 void doCaloJetCorr(edm::Handle<reco::CaloJetCollection> caloJets, edm::Handle<reco::JetCorrector> caloJetCorr);
0076 void doCaloMet(edm::Handle<reco::CaloMETCollection> caloMet);
0077 void doCaloMetBE(edm::Handle<reco::CaloMETCollection> caloMetBE);
0078
0079 void doPFMet(edm::Handle<reco::PFMETCollection> pfMet);
0080 void doPFMetNoMu(edm::Handle<reco::PFMETCollection> pfMet, edm::Handle<reco::MuonCollection>);
0081 void doPUPPIMetNoMu(edm::Handle<reco::PFMETCollection> puppiMet, edm::Handle<reco::MuonCollection>);
0082
0083 void doZPt(edm::Handle<reco::MuonCollection> muons, edm::Handle<std::vector<pat::Jet> > corrPuppiJets);
0084
0085 bool pfJetID(const reco::PFJet& jet);
0086 bool puppiJetID(const pat::Jet& jet);
0087 bool caloJetID(const reco::CaloJet& jet);
0088
0089 public:
0090 L1Analysis::L1AnalysisRecoJetDataFormat* jet_data;
0091 L1Analysis::L1AnalysisRecoMetDataFormat* met_data;
0092
0093 private:
0094
0095 TTree* tree_;
0096
0097
0098 edm::EDGetTokenT<reco::PFJetCollection> pfJetToken_;
0099 edm::EDGetTokenT<reco::PFJetCollection> puppiJetToken_;
0100 edm::EDGetTokenT<std::vector<pat::Jet> > corrPuppiJetToken_;
0101 edm::EDGetTokenT<reco::CaloJetCollection> caloJetToken_;
0102 edm::EDGetTokenT<edm::ValueMap<reco::JetID> > caloJetIDToken_;
0103 edm::EDGetTokenT<reco::JetCorrector> pfJECToken_;
0104 edm::EDGetTokenT<reco::JetCorrector> caloJECToken_;
0105
0106 edm::EDGetTokenT<reco::PFMETCollection> pfMetToken_;
0107 edm::EDGetTokenT<reco::PFMETCollection> puppiMetToken_;
0108 edm::EDGetTokenT<reco::CaloMETCollection> caloMetToken_;
0109 edm::EDGetTokenT<reco::CaloMETCollection> caloMetBEToken_;
0110
0111 edm::EDGetTokenT<reco::MuonCollection> muonToken_;
0112
0113
0114 bool pfJetsMissing_;
0115 bool puppiJetsMissing_;
0116 bool corrPuppiJetsMissing_;
0117 double jetptThreshold_;
0118 double jetetaMax_;
0119 unsigned int maxCl_;
0120 unsigned int maxJet_;
0121 unsigned int maxVtx_;
0122 unsigned int maxTrk_;
0123
0124 bool pfJetCorrMissing_;
0125 bool caloJetCorrMissing_;
0126 bool caloJetsMissing_;
0127 bool caloJetIDMissing_;
0128 bool pfMetMissing_;
0129 bool puppiMetMissing_;
0130 bool caloMetMissing_;
0131 bool caloMetBEMissing_;
0132 bool muonsMissing_;
0133 };
0134
0135 L1JetRecoTreeProducer::L1JetRecoTreeProducer(const edm::ParameterSet& iConfig)
0136 : pfJetsMissing_(false),
0137 puppiJetsMissing_(false),
0138 corrPuppiJetsMissing_(false),
0139 pfJetCorrMissing_(false),
0140 caloJetCorrMissing_(false),
0141 caloJetsMissing_(false),
0142 caloJetIDMissing_(false),
0143 pfMetMissing_(false),
0144 puppiMetMissing_(false),
0145 caloMetMissing_(false),
0146 caloMetBEMissing_(false),
0147 muonsMissing_(false) {
0148 caloJetToken_ =
0149 consumes<reco::CaloJetCollection>(iConfig.getUntrackedParameter("caloJetToken", edm::InputTag("ak4CaloJets")));
0150 pfJetToken_ =
0151 consumes<reco::PFJetCollection>(iConfig.getUntrackedParameter("pfJetToken", edm::InputTag("ak4PFJetsCHS")));
0152 puppiJetToken_ = consumes<reco::PFJetCollection>(iConfig.getParameter<edm::InputTag>("puppiJetToken"));
0153 corrPuppiJetToken_ = consumes<std::vector<pat::Jet> >(
0154 iConfig.getUntrackedParameter("corrPuppiJetToken", edm::InputTag("patJetsCorrectedPuppiJets")));
0155 caloJetIDToken_ =
0156 consumes<edm::ValueMap<reco::JetID> >(iConfig.getUntrackedParameter("caloJetIDToken", edm::InputTag("ak4JetID")));
0157 pfJECToken_ = consumes<reco::JetCorrector>(
0158 iConfig.getUntrackedParameter<edm::InputTag>("pfJECToken", edm::InputTag("ak4PFCHSL1FastL2L3ResidualCorrector")));
0159 caloJECToken_ = consumes<reco::JetCorrector>(iConfig.getUntrackedParameter<edm::InputTag>(
0160 "caloJECToken", edm::InputTag("ak4CaloL1FastL2L3ResidualCorrector")));
0161
0162 pfMetToken_ = consumes<reco::PFMETCollection>(iConfig.getUntrackedParameter("pfMetToken", edm::InputTag("pfMetT1")));
0163 puppiMetToken_ =
0164 consumes<reco::PFMETCollection>(iConfig.getUntrackedParameter("puppiMetToken", edm::InputTag("pfMetPuppi")));
0165 caloMetToken_ =
0166 consumes<reco::CaloMETCollection>(iConfig.getUntrackedParameter("caloMetToken", edm::InputTag("caloMet")));
0167 caloMetBEToken_ =
0168 consumes<reco::CaloMETCollection>(iConfig.getUntrackedParameter("caloMetBEToken", edm::InputTag("caloMetBE")));
0169
0170 muonToken_ = consumes<reco::MuonCollection>(iConfig.getUntrackedParameter("muonToken", edm::InputTag("muons")));
0171
0172 usesResource(TFileService::kSharedResource);
0173
0174 jetptThreshold_ = iConfig.getParameter<double>("jetptThreshold");
0175 jetetaMax_ = iConfig.getParameter<double>("jetetaMax");
0176 maxJet_ = iConfig.getParameter<unsigned int>("maxJet");
0177
0178 jet_data = new L1Analysis::L1AnalysisRecoJetDataFormat();
0179 met_data = new L1Analysis::L1AnalysisRecoMetDataFormat();
0180
0181
0182 edm::Service<TFileService> fs_;
0183 tree_ = fs_->make<TTree>("JetRecoTree", "JetRecoTree");
0184 tree_->Branch("Jet", "L1Analysis::L1AnalysisRecoJetDataFormat", &jet_data, 32000, 3);
0185 tree_->Branch("Sums", "L1Analysis::L1AnalysisRecoMetDataFormat", &met_data, 32000, 3);
0186 }
0187
0188 L1JetRecoTreeProducer::~L1JetRecoTreeProducer() {
0189
0190
0191 }
0192
0193
0194
0195
0196
0197
0198 void L1JetRecoTreeProducer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0199 jet_data->Reset();
0200 met_data->Reset();
0201
0202
0203 edm::Handle<reco::PFJetCollection> pfJets;
0204 iEvent.getByToken(pfJetToken_, pfJets);
0205
0206
0207 edm::Handle<reco::PFJetCollection> puppiJets;
0208 iEvent.getByToken(puppiJetToken_, puppiJets);
0209
0210
0211 edm::Handle<std::vector<pat::Jet> > corrPuppiJets;
0212 iEvent.getByToken(corrPuppiJetToken_, corrPuppiJets);
0213
0214
0215 edm::Handle<reco::CaloJetCollection> caloJets;
0216 iEvent.getByToken(caloJetToken_, caloJets);
0217
0218
0219 edm::Handle<reco::PFMETCollection> pfMet;
0220 iEvent.getByToken(pfMetToken_, pfMet);
0221
0222
0223 edm::Handle<reco::PFMETCollection> puppiMet;
0224 iEvent.getByToken(puppiMetToken_, puppiMet);
0225
0226
0227 edm::Handle<edm::ValueMap<reco::JetID> > jetsID;
0228 iEvent.getByToken(caloJetIDToken_, jetsID);
0229
0230 edm::Handle<reco::JetCorrector> pfJetCorr;
0231 iEvent.getByToken(pfJECToken_, pfJetCorr);
0232
0233 edm::Handle<reco::JetCorrector> caloJetCorr;
0234 iEvent.getByToken(caloJECToken_, caloJetCorr);
0235
0236 edm::Handle<reco::CaloMETCollection> caloMet;
0237 iEvent.getByToken(caloMetToken_, caloMet);
0238
0239 edm::Handle<reco::CaloMETCollection> caloMetBE;
0240 iEvent.getByToken(caloMetBEToken_, caloMetBE);
0241
0242
0243 edm::Handle<reco::MuonCollection> muons;
0244 iEvent.getByToken(muonToken_, muons);
0245
0246 if (pfJets.isValid()) {
0247 jet_data->nJets = 0;
0248
0249 doPFJets(pfJets);
0250
0251 } else {
0252 if (!pfJetsMissing_) {
0253 edm::LogWarning("MissingProduct") << "PFJets not found. Branch will not be filled" << std::endl;
0254 }
0255 pfJetsMissing_ = true;
0256 }
0257
0258 if (pfJetCorr.isValid()) {
0259 doPFJetCorr(pfJets, pfJetCorr);
0260
0261 } else {
0262 if (!pfJetCorrMissing_) {
0263 edm::LogWarning("MissingProduct") << "PF Jet Corrector not found. Branch will not be filled" << std::endl;
0264 }
0265 pfJetCorrMissing_ = true;
0266 }
0267
0268 if (puppiJets.isValid()) {
0269 jet_data->puppi_nUncorrJets = 0;
0270
0271 doPUPPIJets(puppiJets);
0272
0273 } else {
0274 if (!puppiJetsMissing_) {
0275 edm::LogWarning("MissingProduct") << "PUPPIJets not found. Branch will not be filled" << std::endl;
0276 }
0277 puppiJetsMissing_ = true;
0278 }
0279
0280 if (corrPuppiJets.isValid()) {
0281 jet_data->puppi_nJets = 0;
0282
0283 doCorrPUPPIJets(corrPuppiJets);
0284
0285 } else {
0286 if (!corrPuppiJetsMissing_) {
0287 edm::LogWarning("MissingProduct") << "Corrected PUPPIJets not found. Branch will not be filled" << std::endl;
0288 }
0289 corrPuppiJetsMissing_ = true;
0290 }
0291
0292 if (caloJets.isValid()) {
0293 jet_data->nCaloJets = 0;
0294
0295 doCaloJets(caloJets);
0296
0297 } else {
0298 if (!caloJetsMissing_) {
0299 edm::LogWarning("MissingProduct") << "Calo Jets not found. Branch will not be filled" << std::endl;
0300 }
0301 caloJetsMissing_ = true;
0302 }
0303
0304 if (caloJetCorr.isValid()) {
0305 doCaloJetCorr(caloJets, caloJetCorr);
0306
0307 } else {
0308 if (!caloJetCorrMissing_) {
0309 edm::LogWarning("MissingProduct") << "Calo Jet Corrector not found. Branch will not be filled" << std::endl;
0310 }
0311 caloJetCorrMissing_ = true;
0312 }
0313
0314 if (!jetsID.isValid()) {
0315 if (!caloJetIDMissing_) {
0316 edm::LogWarning("MissingProduct") << "Calo Jet ID not found. Branch will not be filled" << std::endl;
0317 }
0318 caloJetIDMissing_ = true;
0319 }
0320
0321 if (pfMet.isValid()) {
0322 doPFMet(pfMet);
0323
0324 if (muons.isValid()) {
0325 doPFMetNoMu(pfMet, muons);
0326
0327 } else {
0328 if (!muonsMissing_) {
0329 edm::LogWarning("MissingProduct") << "Muons not found. PFMetNoMu branch will not be filled" << std::endl;
0330 }
0331 muonsMissing_ = true;
0332 }
0333 } else {
0334 if (!pfMetMissing_) {
0335 edm::LogWarning("MissingProduct") << "PFMet not found. Branch will not be filled" << std::endl;
0336 }
0337 pfMetMissing_ = true;
0338 }
0339
0340 if (puppiMet.isValid()) {
0341 if (muons.isValid()) {
0342 doPUPPIMetNoMu(puppiMet, muons);
0343
0344 } else {
0345 if (!muonsMissing_) {
0346 edm::LogWarning("MissingProduct") << "Muons not found. PUPPIMetNoMu branch will not be filled" << std::endl;
0347 }
0348 muonsMissing_ = true;
0349 }
0350 } else {
0351 if (!puppiMetMissing_) {
0352 edm::LogWarning("MissingProduct") << "PUPPIMet not found. Branch will not be filled" << std::endl;
0353 }
0354 puppiMetMissing_ = true;
0355 }
0356
0357 if (caloMet.isValid()) {
0358 doCaloMet(caloMet);
0359
0360 } else {
0361 if (!caloMetMissing_) {
0362 edm::LogWarning("MissingProduct") << "CaloMet not found. Branch will not be filled" << std::endl;
0363 }
0364 caloMetMissing_ = true;
0365 }
0366
0367 if (caloMetBE.isValid()) {
0368 doCaloMetBE(caloMetBE);
0369
0370 } else {
0371 if (!caloMetBEMissing_) {
0372 edm::LogWarning("MissingProduct") << "CaloMetBE not found. Branch will not be filled" << std::endl;
0373 }
0374 caloMetBEMissing_ = true;
0375 }
0376
0377 if (muons.isValid()) {
0378 if (puppiJets.isValid()) {
0379 doZPt(muons, corrPuppiJets);
0380
0381 } else {
0382 if (!puppiJetsMissing_) {
0383 edm::LogWarning("MissingProduct") << "PUPPIJets not found. Branch will not be filled" << std::endl;
0384 }
0385 puppiJetsMissing_ = true;
0386 }
0387 } else {
0388 if (!muonsMissing_) {
0389 edm::LogWarning("MissingProduct") << "Muons not found. ZPt branch will not be filled" << std::endl;
0390 }
0391 muonsMissing_ = true;
0392 }
0393
0394 tree_->Fill();
0395 }
0396
0397 void L1JetRecoTreeProducer::doCaloJets(edm::Handle<reco::CaloJetCollection> caloJets) {
0398 for (auto it = caloJets->begin(); it != caloJets->end() && jet_data->nCaloJets < maxJet_; ++it) {
0399 if (!caloJetIDMissing_)
0400 if (!caloJetID(*it))
0401 continue;
0402
0403 jet_data->caloEt.push_back(it->et());
0404 jet_data->caloEta.push_back(it->eta());
0405 jet_data->caloPhi.push_back(it->phi());
0406 jet_data->caloE.push_back(it->energy());
0407
0408 jet_data->eEMF.push_back(it->emEnergyFraction());
0409 jet_data->eEmEB.push_back(it->emEnergyInEB());
0410 jet_data->eEmEE.push_back(it->emEnergyInEE());
0411 jet_data->eEmHF.push_back(it->emEnergyInHF());
0412 jet_data->eHadHB.push_back(it->hadEnergyInHB());
0413 jet_data->eHadHE.push_back(it->hadEnergyInHE());
0414 jet_data->eHadHO.push_back(it->hadEnergyInHO());
0415 jet_data->eHadHF.push_back(it->hadEnergyInHF());
0416 jet_data->eMaxEcalTow.push_back(it->maxEInEmTowers());
0417 jet_data->eMaxHcalTow.push_back(it->maxEInHadTowers());
0418 jet_data->towerArea.push_back(it->towersArea());
0419 jet_data->n60.push_back(it->n60());
0420
0421 jet_data->nCaloJets++;
0422 }
0423 }
0424
0425 void L1JetRecoTreeProducer::doPFJets(edm::Handle<reco::PFJetCollection> pfJets) {
0426 for (auto it = pfJets->begin(); it != pfJets->end() && jet_data->nJets < maxJet_; ++it) {
0427 if (!pfJetID(*it))
0428 continue;
0429
0430 jet_data->et.push_back(it->et());
0431 jet_data->eta.push_back(it->eta());
0432 jet_data->phi.push_back(it->phi());
0433 jet_data->e.push_back(it->energy());
0434
0435 jet_data->chef.push_back(it->chargedHadronEnergyFraction());
0436 jet_data->nhef.push_back(it->neutralHadronEnergyFraction());
0437 jet_data->pef.push_back(it->photonEnergyFraction());
0438 jet_data->eef.push_back(it->electronEnergyFraction());
0439 jet_data->mef.push_back(it->muonEnergyFraction());
0440 jet_data->hfhef.push_back(it->HFHadronEnergyFraction());
0441 jet_data->hfemef.push_back(it->HFEMEnergyFraction());
0442 jet_data->chMult.push_back(it->chargedHadronMultiplicity());
0443 jet_data->nhMult.push_back(it->neutralHadronMultiplicity());
0444 jet_data->phMult.push_back(it->photonMultiplicity());
0445 jet_data->elMult.push_back(it->electronMultiplicity());
0446 jet_data->muMult.push_back(it->muonMultiplicity());
0447 jet_data->hfhMult.push_back(it->HFHadronMultiplicity());
0448 jet_data->hfemMult.push_back(it->HFEMMultiplicity());
0449
0450 jet_data->cemef.push_back(it->chargedEmEnergyFraction());
0451 jet_data->cmef.push_back(it->chargedMuEnergyFraction());
0452 jet_data->nemef.push_back(it->neutralEmEnergyFraction());
0453 jet_data->cMult.push_back(it->chargedMultiplicity());
0454 jet_data->nMult.push_back(it->neutralMultiplicity());
0455
0456 jet_data->nJets++;
0457 }
0458 }
0459
0460 void L1JetRecoTreeProducer::doPFJetCorr(edm::Handle<reco::PFJetCollection> pfJets,
0461 edm::Handle<reco::JetCorrector> pfJetCorr) {
0462 float corrFactor = 1.;
0463 unsigned int nJets = 0;
0464
0465 float mHx = 0;
0466 float mHy = 0;
0467
0468 met_data->Ht = 0;
0469 met_data->mHt = -999.;
0470 met_data->mHtPhi = -999.;
0471
0472 for (auto it = pfJets->begin(); it != pfJets->end() && nJets < maxJet_; ++it) {
0473 if (!pfJetID(*it))
0474 continue;
0475
0476 corrFactor = pfJetCorr.product()->correction(*it);
0477
0478 jet_data->etCorr.push_back(it->et() * corrFactor);
0479 jet_data->corrFactor.push_back(corrFactor);
0480
0481 nJets++;
0482
0483 if (it->pt() * corrFactor > jetptThreshold_ && std::abs(it->eta()) < jetetaMax_) {
0484 mHx += -1. * it->px() * corrFactor;
0485 mHy += -1. * it->py() * corrFactor;
0486 met_data->Ht += it->pt() * corrFactor;
0487 }
0488 }
0489
0490 TVector2 tv2 = TVector2(mHx, mHy);
0491 met_data->mHt = tv2.Mod();
0492 met_data->mHtPhi = tv2.Phi();
0493 }
0494
0495 void L1JetRecoTreeProducer::doPUPPIJets(edm::Handle<reco::PFJetCollection> puppiJets) {
0496 for (auto it = puppiJets->begin(); it != puppiJets->end() && jet_data->puppi_nUncorrJets < maxJet_; ++it) {
0497 if (!puppiJetID(*it))
0498 continue;
0499 jet_data->puppi_et.push_back(it->et());
0500 jet_data->puppi_nUncorrJets++;
0501 }
0502 }
0503
0504 void L1JetRecoTreeProducer::doCorrPUPPIJets(edm::Handle<std::vector<pat::Jet> > corrPuppiJets) {
0505 float mHx = 0;
0506 float mHy = 0;
0507
0508 met_data->puppi_Ht = 0;
0509 met_data->puppi_mHt = -999.;
0510 met_data->puppi_mHtPhi = -999.;
0511
0512 for (auto it = corrPuppiJets->begin(); it != corrPuppiJets->end() && jet_data->puppi_nJets < maxJet_; ++it) {
0513 if (!puppiJetID(*it))
0514 continue;
0515
0516 jet_data->puppi_etCorr.push_back(it->et());
0517 jet_data->puppi_eta.push_back(it->eta());
0518 jet_data->puppi_phi.push_back(it->phi());
0519 jet_data->puppi_e.push_back(it->energy());
0520
0521 jet_data->puppi_chef.push_back(it->chargedHadronEnergyFraction());
0522 jet_data->puppi_nhef.push_back(it->neutralHadronEnergyFraction());
0523 jet_data->puppi_pef.push_back(it->photonEnergyFraction());
0524 jet_data->puppi_eef.push_back(it->electronEnergyFraction());
0525 jet_data->puppi_mef.push_back(it->muonEnergyFraction());
0526 jet_data->puppi_hfhef.push_back(it->HFHadronEnergyFraction());
0527 jet_data->puppi_hfemef.push_back(it->HFEMEnergyFraction());
0528 jet_data->puppi_chMult.push_back(it->chargedHadronMultiplicity());
0529 jet_data->puppi_nhMult.push_back(it->neutralHadronMultiplicity());
0530 jet_data->puppi_phMult.push_back(it->photonMultiplicity());
0531 jet_data->puppi_elMult.push_back(it->electronMultiplicity());
0532 jet_data->puppi_muMult.push_back(it->muonMultiplicity());
0533 jet_data->puppi_hfhMult.push_back(it->HFHadronMultiplicity());
0534 jet_data->puppi_hfemMult.push_back(it->HFEMMultiplicity());
0535
0536 jet_data->puppi_cemef.push_back(it->chargedEmEnergyFraction());
0537 jet_data->puppi_cmef.push_back(it->chargedMuEnergyFraction());
0538 jet_data->puppi_nemef.push_back(it->neutralEmEnergyFraction());
0539 jet_data->puppi_cMult.push_back(it->chargedMultiplicity());
0540 jet_data->puppi_nMult.push_back(it->neutralMultiplicity());
0541
0542 jet_data->puppi_nJets++;
0543
0544 if (it->pt() > jetptThreshold_ && std::abs(it->eta()) < jetetaMax_) {
0545 mHx += -1. * it->px();
0546 mHy += -1. * it->py();
0547 met_data->puppi_Ht += it->pt();
0548 }
0549 }
0550
0551 TVector2 tv2 = TVector2(mHx, mHy);
0552 met_data->puppi_mHt = tv2.Mod();
0553 met_data->puppi_mHtPhi = tv2.Phi();
0554 }
0555
0556 void L1JetRecoTreeProducer::doCaloJetCorr(edm::Handle<reco::CaloJetCollection> caloJets,
0557 edm::Handle<reco::JetCorrector> caloJetCorr) {
0558 float caloCorrFactor = 1.;
0559 unsigned int nCaloJets = 0;
0560
0561 met_data->caloHt = 0;
0562
0563 for (auto it = caloJets->begin(); it != caloJets->end() && nCaloJets < maxJet_; ++it) {
0564 if (!caloJetIDMissing_)
0565 if (!caloJetID(*it))
0566 continue;
0567
0568 caloCorrFactor = caloJetCorr.product()->correction(*it);
0569
0570 jet_data->caloEtCorr.push_back(it->et() * caloCorrFactor);
0571 jet_data->caloCorrFactor.push_back(caloCorrFactor);
0572
0573 nCaloJets++;
0574
0575 if (it->pt() * caloCorrFactor > jetptThreshold_ && std::abs(it->eta()) < jetetaMax_) {
0576 met_data->caloHt += it->pt() * caloCorrFactor;
0577 }
0578 }
0579 }
0580
0581 void L1JetRecoTreeProducer::doPFMet(edm::Handle<reco::PFMETCollection> pfMet) {
0582 const reco::PFMETCollection* metCol = pfMet.product();
0583 const reco::PFMET theMet = metCol->front();
0584
0585 met_data->met = theMet.et();
0586 met_data->metPhi = theMet.phi();
0587 met_data->sumEt = theMet.sumEt();
0588 met_data->metPx = theMet.px();
0589 met_data->metPy = theMet.py();
0590 }
0591
0592 void L1JetRecoTreeProducer::doPFMetNoMu(edm::Handle<reco::PFMETCollection> pfMet,
0593 edm::Handle<reco::MuonCollection> muons) {
0594 const reco::PFMETCollection* metCol = pfMet.product();
0595 const reco::PFMET theMet = metCol->front();
0596 reco::PFMET thePFMetNoMu = metCol->front();
0597
0598 double pfMetNoMuPx = theMet.px();
0599 double pfMetNoMuPy = theMet.py();
0600
0601 double muPx(0.), muPy(0.);
0602
0603 for (auto it = muons->begin(); it != muons->end(); ++it) {
0604 if (it->isPFMuon()) {
0605 muPx += it->px();
0606 muPy += it->py();
0607 }
0608 }
0609
0610 pfMetNoMuPx += muPx;
0611 pfMetNoMuPy += muPy;
0612
0613 math::XYZTLorentzVector pfMetNoMuP4(pfMetNoMuPx, pfMetNoMuPy, 0, hypot(pfMetNoMuPx, pfMetNoMuPy));
0614
0615 thePFMetNoMu.setP4(pfMetNoMuP4);
0616
0617 met_data->pfMetNoMu = thePFMetNoMu.et();
0618 met_data->pfMetNoMuPhi = thePFMetNoMu.phi();
0619 met_data->pfMetNoMuPx = thePFMetNoMu.px();
0620 met_data->pfMetNoMuPy = thePFMetNoMu.py();
0621 }
0622
0623 void L1JetRecoTreeProducer::doPUPPIMetNoMu(edm::Handle<reco::PFMETCollection> puppiMet,
0624 edm::Handle<reco::MuonCollection> muons) {
0625 const reco::PFMETCollection* metCol = puppiMet.product();
0626 const reco::PFMET theMet = metCol->front();
0627 reco::PFMET thePUPPIMetNoMu = metCol->front();
0628
0629 double puppiMetNoMuPx = theMet.px();
0630 double puppiMetNoMuPy = theMet.py();
0631
0632 double muPx(0.), muPy(0.);
0633
0634 for (auto it = muons->begin(); it != muons->end(); ++it) {
0635 if (it->isPFMuon()) {
0636 muPx += it->px();
0637 muPy += it->py();
0638 }
0639 }
0640
0641 puppiMetNoMuPx += muPx;
0642 puppiMetNoMuPy += muPy;
0643
0644 math::XYZTLorentzVector puppiMetNoMuP4(puppiMetNoMuPx, puppiMetNoMuPy, 0, hypot(puppiMetNoMuPx, puppiMetNoMuPy));
0645
0646 thePUPPIMetNoMu.setP4(puppiMetNoMuP4);
0647
0648 met_data->puppi_metNoMu = thePUPPIMetNoMu.et();
0649 met_data->puppi_metNoMuPhi = thePUPPIMetNoMu.phi();
0650 met_data->puppi_metNoMuPx = thePUPPIMetNoMu.px();
0651 met_data->puppi_metNoMuPy = thePUPPIMetNoMu.py();
0652 }
0653
0654 void L1JetRecoTreeProducer::doCaloMet(edm::Handle<reco::CaloMETCollection> caloMet) {
0655 const reco::CaloMETCollection* metCol = caloMet.product();
0656 const reco::CaloMET theMet = metCol->front();
0657
0658 met_data->caloMet = theMet.et();
0659 met_data->caloMetPhi = theMet.phi();
0660 met_data->caloSumEt = theMet.sumEt();
0661 }
0662
0663 void L1JetRecoTreeProducer::doCaloMetBE(edm::Handle<reco::CaloMETCollection> caloMetBE) {
0664 const reco::CaloMETCollection* metCol = caloMetBE.product();
0665 const reco::CaloMET theMet = metCol->front();
0666
0667 met_data->caloMetBE = theMet.et();
0668 met_data->caloMetPhiBE = theMet.phi();
0669 met_data->caloSumEtBE = theMet.sumEt();
0670 }
0671
0672 void L1JetRecoTreeProducer::doZPt(edm::Handle<reco::MuonCollection> muons,
0673 edm::Handle<std::vector<pat::Jet> > corrPuppiJets) {
0674 bool passPuppiJetPtCut = false;
0675
0676 for (auto it = corrPuppiJets->begin(); it != corrPuppiJets->end(); ++it) {
0677 if (!puppiJetID(*it))
0678 continue;
0679 if (it->muonEnergyFraction() > 0.5 || it->chargedEmEnergyFraction() > 0.5)
0680 continue;
0681 if (it->pt() > 30)
0682 passPuppiJetPtCut = true;
0683 }
0684
0685 if (!passPuppiJetPtCut || muons->size() < 2) {
0686 met_data->zPt = -999;
0687 return;
0688 }
0689
0690 reco::Muon muon1;
0691 reco::Muon muon2;
0692
0693 float zMass = 91.2;
0694 float diMuMass = 0;
0695 float closestDiff = 999.;
0696 bool found2PFMuons = false;
0697
0698 for (auto it1 = muons->begin(); it1 != muons->end(); ++it1) {
0699 if (!it1->isPFMuon())
0700 continue;
0701 for (auto it2 = std::next(it1); it2 != muons->end(); ++it2) {
0702 if (!it2->isPFMuon())
0703 continue;
0704 if (it1->charge() != (-1 * it2->charge()))
0705 continue;
0706
0707 found2PFMuons = true;
0708 diMuMass = (it1->p4() + it2->p4()).M();
0709 float diff = abs(diMuMass - zMass);
0710 if (diff < closestDiff) {
0711 closestDiff = diff;
0712 muon1 = *it1;
0713 muon2 = *it2;
0714 }
0715 }
0716 }
0717
0718 diMuMass = (muon1.p4() + muon2.p4()).M();
0719 if (abs(diMuMass - zMass) > 30 || !found2PFMuons) {
0720 met_data->zPt = -999;
0721 return;
0722 }
0723
0724 float zPt = (muon1.p4() + muon2.p4()).pt();
0725 met_data->zPt = zPt;
0726 }
0727
0728 bool L1JetRecoTreeProducer::pfJetID(const reco::PFJet& jet) {
0729 bool tmp = true;
0730 if (std::abs(jet.eta()) <= 2.6) {
0731 tmp &= jet.neutralHadronEnergyFraction() < 0.9;
0732 tmp &= jet.neutralEmEnergyFraction() < 0.9;
0733 tmp &= (jet.chargedMultiplicity() + jet.neutralMultiplicity()) > 1;
0734 tmp &= jet.muonEnergyFraction() < 0.8;
0735 tmp &= jet.chargedHadronEnergyFraction() > 0.01;
0736 tmp &= jet.chargedMultiplicity() > 0;
0737 tmp &= jet.chargedEmEnergyFraction() < 0.8;
0738 }
0739 if (std::abs(jet.eta()) > 2.6 && std::abs(jet.eta()) <= 2.7) {
0740 tmp &= jet.neutralHadronEnergyFraction() < 0.9;
0741 tmp &= jet.neutralEmEnergyFraction() < 0.99;
0742 tmp &= jet.muonEnergyFraction() < 0.8;
0743 tmp &= jet.chargedMultiplicity() > 0;
0744 tmp &= jet.chargedEmEnergyFraction() < 0.8;
0745 }
0746 if (std::abs(jet.eta()) > 2.7 && std::abs(jet.eta()) < 3.0) {
0747 tmp &= jet.neutralEmEnergyFraction() < 0.99;
0748 tmp &= jet.neutralMultiplicity() > 1;
0749 }
0750 if (std::abs(jet.eta()) > 3.0) {
0751 tmp &= jet.neutralHadronEnergyFraction() > 0.2;
0752 tmp &= jet.neutralEmEnergyFraction() < 0.9;
0753 tmp &= jet.neutralMultiplicity() > 10;
0754 }
0755
0756
0757
0758
0759
0760 return tmp;
0761 }
0762
0763
0764 bool L1JetRecoTreeProducer::puppiJetID(const pat::Jet& jet) {
0765 bool tmp = true;
0766 if (std::abs(jet.eta()) <= 2.6) {
0767 tmp &= jet.neutralHadronEnergyFraction() < 0.9;
0768 tmp &= jet.neutralEmEnergyFraction() < 0.9;
0769 tmp &= (jet.chargedMultiplicity() + jet.neutralMultiplicity()) > 1;
0770 tmp &= jet.muonEnergyFraction() < 0.8;
0771 tmp &= jet.chargedHadronEnergyFraction() > 0.01;
0772 tmp &= jet.chargedMultiplicity() > 0;
0773 tmp &= jet.chargedEmEnergyFraction() < 0.8;
0774 }
0775 if (std::abs(jet.eta()) > 2.6 && std::abs(jet.eta()) <= 2.7) {
0776 tmp &= jet.neutralHadronEnergyFraction() < 0.9;
0777 tmp &= jet.neutralEmEnergyFraction() < 0.99;
0778 tmp &= jet.muonEnergyFraction() < 0.8;
0779 tmp &= jet.chargedEmEnergyFraction() < 0.8;
0780 }
0781 if (std::abs(jet.eta()) > 2.7 && std::abs(jet.eta()) <= 3.0) {
0782 tmp &= jet.neutralHadronEnergyFraction() < 0.9999;
0783 }
0784 if (std::abs(jet.eta()) > 3.0) {
0785 tmp &= jet.neutralEmEnergyFraction() < 0.9;
0786 tmp &= jet.neutralMultiplicity() > 2;
0787 }
0788
0789 return tmp;
0790 }
0791
0792 bool L1JetRecoTreeProducer::caloJetID(const reco::CaloJet& jet) {
0793 bool tmp = true;
0794
0795 return tmp;
0796 }
0797
0798
0799 void L1JetRecoTreeProducer::beginJob(void) {}
0800
0801
0802 void L1JetRecoTreeProducer::endJob() {}
0803
0804
0805 DEFINE_FWK_MODULE(L1JetRecoTreeProducer);
