Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-01-24 01:11:28

0001 #include "DQM/Physics/src/ExoticaDQM.h"
0002 
0003 using namespace edm;
0004 using namespace std;
0005 using namespace reco;
0006 using namespace trigger;
0007 
0008 typedef vector<string> vstring;
0009 
0010 //
0011 // -- Constructor
0012 //
0013 ExoticaDQM::ExoticaDQM(const edm::ParameterSet& ps) {
0014   edm::LogInfo("ExoticaDQM") << " Starting ExoticaDQM "
0015                              << "\n";
0016 
0017   typedef std::vector<edm::InputTag> vtag;
0018 
0019   // Get parameters from configuration file
0020   // Trigger
0021   TriggerToken_ = consumes<TriggerResults>(ps.getParameter<edm::InputTag>("TriggerResults"));
0022   HltPaths_ = ps.getParameter<vector<string> >("HltPaths");
0023   //
0024   VertexToken_ = consumes<reco::VertexCollection>(ps.getParameter<InputTag>("vertexCollection"));
0025   //
0026   ElectronToken_ = consumes<reco::GsfElectronCollection>(ps.getParameter<InputTag>("electronCollection"));
0027   //
0028   MuonToken_ = consumes<reco::MuonCollection>(ps.getParameter<InputTag>("muonCollection"));
0029   //
0030   PhotonToken_ = consumes<reco::PhotonCollection>(ps.getParameter<InputTag>("photonCollection"));
0031   //
0032   PFJetToken_ = consumes<reco::PFJetCollection>(ps.getParameter<InputTag>("pfJetCollection"));
0033   //
0034   DiJetPFJetCollection_ = ps.getParameter<std::vector<edm::InputTag> >("DiJetPFJetCollection");
0035   for (std::vector<edm::InputTag>::const_iterator jetlabel = DiJetPFJetCollection_.begin(),
0036                                                   jetlabelEnd = DiJetPFJetCollection_.end();
0037        jetlabel != jetlabelEnd;
0038        ++jetlabel) {
0039     DiJetPFJetToken_.push_back(consumes<reco::PFJetCollection>(*jetlabel));
0040   }
0041   //
0042   PFMETToken_ = consumes<reco::PFMETCollection>(ps.getParameter<InputTag>("pfMETCollection"));
0043   //
0044   ecalBarrelRecHitToken_ = consumes<EBRecHitCollection>(
0045       ps.getUntrackedParameter<InputTag>("ecalBarrelRecHit", InputTag("reducedEcalRecHitsEB")));
0046   //
0047   ecalEndcapRecHitToken_ = consumes<EERecHitCollection>(
0048       ps.getUntrackedParameter<InputTag>("ecalEndcapRecHit", InputTag("reducedEcalRecHitsEE")));
0049   //
0050   TrackToken_ = consumes<reco::TrackCollection>(ps.getParameter<InputTag>("trackCollection"));
0051   //
0052   MuonDispToken_ = consumes<reco::TrackCollection>(ps.getParameter<InputTag>("displacedMuonCollection"));
0053   //
0054   MuonDispSAToken_ = consumes<reco::TrackCollection>(ps.getParameter<InputTag>("displacedSAMuonCollection"));
0055   //
0056   GenParticleToken_ = consumes<reco::GenParticleCollection>(ps.getParameter<InputTag>("genParticleCollection"));
0057 
0058   JetCorrectorToken_ = consumes<reco::JetCorrector>(ps.getParameter<edm::InputTag>("jetCorrector"));
0059 
0060   magFieldToken_ = esConsumes();
0061   //Cuts - MultiJets
0062   jetID = new reco::helper::JetIDHelper(ps.getParameter<ParameterSet>("JetIDParams"), consumesCollector());
0063 
0064   //Varibles and Cuts for each Module:
0065   //Dijet
0066   dijet_PFJet1_pt_cut_ = ps.getParameter<double>("dijet_PFJet1_pt_cut");
0067   dijet_PFJet2_pt_cut_ = ps.getParameter<double>("dijet_PFJet2_pt_cut");
0068   //DiMuon
0069   dimuon_Muon1_pt_cut_ = ps.getParameter<double>("dimuon_Muon1_pt_cut");
0070   dimuon_Muon2_pt_cut_ = ps.getParameter<double>("dimuon_Muon2_pt_cut");
0071   //DiElectron
0072   dielectron_Electron1_pt_cut_ = ps.getParameter<double>("dielectron_Electron2_pt_cut");
0073   dielectron_Electron2_pt_cut_ = ps.getParameter<double>("dielectron_Electron2_pt_cut");
0074   //DiPhoton
0075   diphoton_Photon1_pt_cut_ = ps.getParameter<double>("diphoton_Photon2_pt_cut");
0076   diphoton_Photon2_pt_cut_ = ps.getParameter<double>("diphoton_Photon2_pt_cut");
0077   //MonoJet
0078   monojet_PFJet_pt_cut_ = ps.getParameter<double>("monojet_PFJet_pt_cut");
0079   monojet_PFJet_met_cut_ = ps.getParameter<double>("monojet_PFJet_met_cut");
0080   //MonoMuon
0081   monomuon_Muon_pt_cut_ = ps.getParameter<double>("monomuon_Muon_pt_cut");
0082   monomuon_Muon_met_cut_ = ps.getParameter<double>("monomuon_Muon_met_cut");
0083   //MonoElectron
0084   monoelectron_Electron_pt_cut_ = ps.getParameter<double>("monoelectron_Electron_pt_cut");
0085   monoelectron_Electron_met_cut_ = ps.getParameter<double>("monoelectron_Electron_met_cut");
0086   //MonoPhoton
0087   monophoton_Photon_pt_cut_ = ps.getParameter<double>("monophoton_Photon_pt_cut");
0088   monophoton_Photon_met_cut_ = ps.getParameter<double>("monophoton_Photon_met_cut");
0089   // Displaced lepton or jet
0090   dispFermion_eta_cut_ = ps.getParameter<double>("dispFermion_eta_cut");
0091   dispFermion_pt_cut_ = ps.getParameter<double>("dispFermion_pt_cut");
0092 }
0093 
0094 //
0095 // -- Destructor
0096 //
0097 ExoticaDQM::~ExoticaDQM() {
0098   edm::LogInfo("ExoticaDQM") << " Deleting ExoticaDQM "
0099                              << "\n";
0100 }
0101 
0102 //
0103 //  -- Book histograms
0104 //
0105 void ExoticaDQM::bookHistograms(DQMStore::IBooker& bei, edm::Run const&, edm::EventSetup const&) {
0106   bei.cd();
0107 
0108   //--- DiJet
0109   for (unsigned int icoll = 0; icoll < DiJetPFJetCollection_.size(); ++icoll) {
0110     std::stringstream ss;
0111     ss << "Physics/Exotica/Dijets/" << DiJetPFJetCollection_[icoll].label();
0112     bei.setCurrentFolder(ss.str());
0113     //bei.setCurrentFolder("Physics/Exotica/Dijets");
0114     dijet_PFJet_pt.push_back(bei.book1D("dijet_PFJet_pt", "Pt of PFJet (GeV)", 50, 30.0, 5000));
0115     dijet_PFJet_eta.push_back(bei.book1D("dijet_PFJet_eta", "#eta (PFJet)", 50, -2.5, 2.5));
0116     dijet_PFJet_phi.push_back(bei.book1D("dijet_PFJet_phi", "#phi (PFJet)", 50, -3.14, 3.14));
0117     dijet_PFJet_rapidity.push_back(bei.book1D("dijet_PFJet_rapidity", "Rapidity (PFJet)", 50, -6.0, 6.0));
0118     dijet_PFJet_mass.push_back(bei.book1D("dijet_PFJet_mass", "Mass (PFJets)", 50, 0., 500.));
0119     dijet_deltaPhiPFJet1PFJet2.push_back(
0120         bei.book1D("dijet_deltaPhiPFJet1PFJet2", "#Delta#phi(Leading PFJet, Sub PFJet)", 40, 0., 3.15));
0121     dijet_deltaEtaPFJet1PFJet2.push_back(
0122         bei.book1D("dijet_deltaEtaPFJet1PFJet2", "#Delta#eta(Leading PFJet, Sub PFJet)", 40, -5., 5.));
0123     dijet_deltaRPFJet1PFJet2.push_back(
0124         bei.book1D("dijet_deltaRPFJet1PFJet2", "#DeltaR(Leading PFJet, Sub PFJet)", 50, 0., 6.));
0125     dijet_invMassPFJet1PFJet2.push_back(
0126         bei.book1D("dijet_invMassPFJet1PFJet2", "Leading PFJet, SubLeading PFJet Invariant mass (GeV)", 50, 0., 8000.));
0127     dijet_PFchef.push_back(bei.book1D("dijet_PFchef", "Leading PFJet CHEF", 50, 0.0, 1.0));
0128     dijet_PFnhef.push_back(bei.book1D("dijet_PFnhef", "Leading PFJet NHEF", 50, 0.0, 1.0));
0129     dijet_PFcemf.push_back(bei.book1D("dijet_PFcemf", "Leading PFJet CEMF", 50, 0.0, 1.0));
0130     dijet_PFnemf.push_back(bei.book1D("dijet_PFnemf", "Leading PFJEt NEMF", 50, 0.0, 1.0));
0131     dijet_PFJetMulti.push_back(bei.book1D("dijet_PFJetMulti", "No. of PFJets", 10, 0., 10.));
0132   }
0133   //--- DiMuon
0134   bei.setCurrentFolder("Physics/Exotica/DiMuons");
0135   dimuon_Muon_pt = bei.book1D("dimuon_Muon_pt", "Pt of Muon (GeV)", 50, 30.0, 2000);
0136   dimuon_Muon_eta = bei.book1D("dimuon_Muon_eta", "#eta (Muon)", 50, -2.5, 2.5);
0137   dimuon_Muon_phi = bei.book1D("dimuon_Muon_phi", "#phi (Muon)", 50, -3.14, 3.14);
0138   dimuon_Charge = bei.book1D("dimuon_Charge", "Charge of the Muon", 10, -5., 5.);
0139   dimuon_deltaEtaMuon1Muon2 =
0140       bei.book1D("dimuon_deltaEtaMuon1Muon2", "#Delta#eta(Leading Muon, Sub Muon)", 40, -5., 5.);
0141   dimuon_deltaPhiMuon1Muon2 =
0142       bei.book1D("dimuon_deltaPhiMuon1Muon2", "#Delta#phi(Leading Muon, Sub Muon)", 40, 0., 3.15);
0143   dimuon_deltaRMuon1Muon2 = bei.book1D("dimuon_deltaRMuon1Muon2", "#DeltaR(Leading Muon, Sub Muon)", 50, 0., 6.);
0144   dimuon_invMassMuon1Muon2 =
0145       bei.book1D("dimuon_invMassMuon1Muon2", "Leading Muon, SubLeading Muon Low Invariant mass (GeV)", 50, 500., 4500.);
0146   dimuon_MuonMulti = bei.book1D("dimuon_MuonMulti", "No. of Muons", 10, 0., 10.);
0147   //--- DiElectrons
0148   bei.setCurrentFolder("Physics/Exotica/DiElectrons");
0149   dielectron_Electron_pt = bei.book1D("dielectron_Electron_pt", "Pt of Electron (GeV)", 50, 30.0, 2000);
0150   dielectron_Electron_eta = bei.book1D("dielectron_Electron_eta", "#eta (Electron)", 50, -2.5, 2.5);
0151   dielectron_Electron_phi = bei.book1D("dielectron_Electron_phi", "#phi (Electron)", 50, -3.14, 3.14);
0152   dielectron_Charge = bei.book1D("dielectron_Charge", "Charge of the Electron", 10, -5., 5.);
0153   dielectron_deltaEtaElectron1Electron2 =
0154       bei.book1D("dielectron_deltaEtaElectron1Electron2", "#Delta#eta(Leading Electron, Sub Electron)", 40, -5., 5.);
0155   dielectron_deltaPhiElectron1Electron2 =
0156       bei.book1D("dielectron_deltaPhiElectron1Electron2", "#Delta#phi(Leading Electron, Sub Electron)", 40, 0., 3.15);
0157   dielectron_deltaRElectron1Electron2 =
0158       bei.book1D("dielectron_deltaRElectron1Electron2", "#DeltaR(Leading Electron, Sub Electron)", 50, 0., 6.);
0159   dielectron_invMassElectron1Electron2 = bei.book1D("dielectron_invMassElectron1Electron2",
0160                                                     "Leading Electron, SubLeading Electron Invariant mass (GeV)",
0161                                                     50,
0162                                                     500.,
0163                                                     4500.);
0164   dielectron_ElectronMulti = bei.book1D("dielectron_ElectronMulti", "No. of Electrons", 10, 0., 10.);
0165   //--- DiPhotons
0166   bei.setCurrentFolder("Physics/Exotica/DiPhotons");
0167   diphoton_Photon_energy = bei.book1D("diphoton_Photon_energy", "Energy of Photon (GeV)", 50, 30.0, 300);
0168   diphoton_Photon_et = bei.book1D("diphoton_Photon_et", "Et of Photon (GeV)", 50, 30.0, 300);
0169   diphoton_Photon_pt = bei.book1D("diphoton_Photon_pt", "Pt of Photon (GeV)", 50, 30.0, 300);
0170   diphoton_Photon_eta = bei.book1D("diphoton_Photon_eta", "#eta (Photon)", 50, -2.5, 2.5);
0171   diphoton_Photon_etasc = bei.book1D("diphoton_Photon_etasc", "#eta sc(Photon)", 50, -2.5, 2.5);
0172   diphoton_Photon_phi = bei.book1D("diphoton_Photon_phi", "#phi (Photon)", 50, -3.14, 3.14);
0173   diphoton_Photon_hovere_eb = bei.book1D("diphoton_Photon_hovere_eb", "H/E (Photon) EB", 50, 0., 0.50);
0174   diphoton_Photon_hovere_ee = bei.book1D("diphoton_Photon_hovere_ee", "H/E (Photon) EE", 50, 0., 0.50);
0175   diphoton_Photon_sigmaietaieta_eb =
0176       bei.book1D("diphoton_Photon_sigmaietaieta_eb", "#sigma_{i #eta i #eta} (Photon) EB", 50, 0., 0.03);
0177   diphoton_Photon_sigmaietaieta_ee =
0178       bei.book1D("diphoton_Photon_sigmaietaieta_ee", "#sigma_{i #eta i #eta} (Photon) EE", 50, 0., 0.03);
0179   diphoton_Photon_trksumptsolidconedr03_eb =
0180       bei.book1D("diphoton_Photon_trksumptsolidconedr03_eb", "TrkSumPtDr03 (Photon) EB", 50, 0., 15.);
0181   diphoton_Photon_trksumptsolidconedr03_ee =
0182       bei.book1D("diphoton_Photon_trksumptsolidconedr03_ee", "TrkSumPtDr03 (Photon) EE", 50, 0., 15.);
0183   diphoton_Photon_e1x5e5x5_eb = bei.book1D("diphoton_Photon_e1x5e5x5_eb", "E_{1x5}/E_{5x5} (Photon) EB", 50, 0., 1.);
0184   diphoton_Photon_e1x5e5x5_ee = bei.book1D("diphoton_Photon_e1x5e5x5_ee", "E_{1x5}/E_{5x5} (Photon) EE", 50, 0., 1.);
0185   diphoton_Photon_e2x5e5x5_eb = bei.book1D("diphoton_Photon_e2x5e5x5_eb", "E_{2x5}/E_{5x5} (Photon) EB", 50, 0., 1.);
0186   diphoton_Photon_e2x5e5x5_ee = bei.book1D("diphoton_Photon_e2x5e5x5_ee", "E_{2x5}/E_{5x5} (Photon) EE", 50, 0., 1.);
0187   diphoton_deltaEtaPhoton1Photon2 =
0188       bei.book1D("diphoton_deltaEtaPhoton1Photon2", "#Delta#eta(SubLeading Photon, Sub Photon)", 40, -5., 5.);
0189   diphoton_deltaPhiPhoton1Photon2 =
0190       bei.book1D("diphoton_deltaPhiPhoton1Photon2", "#Delta#phi(SubLeading Photon, Sub Photon)", 40, 0., 3.15);
0191   diphoton_deltaRPhoton1Photon2 =
0192       bei.book1D("diphoton_deltaRPhoton1Photon2", "#DeltaR(SubLeading Photon, Sub Photon)", 50, 0., 6.);
0193   diphoton_invMassPhoton1Photon2 = bei.book1D(
0194       "diphoton_invMassPhoton1Photon2", "SubLeading Photon, SubSubLeading Photon Invariant mass (GeV)", 50, 500., 4500.);
0195   diphoton_PhotonMulti = bei.book1D("diphoton_PhotonMulti", "No. of Photons", 10, 0., 10.);
0196   //--- MonoJet
0197   bei.setCurrentFolder("Physics/Exotica/MonoJet");
0198   monojet_PFJet_pt = bei.book1D("monojet_PFJet_pt", "Pt of MonoJet (GeV)", 50, 30.0, 1000);
0199   monojet_PFJet_eta = bei.book1D("monojet_PFJet_eta", "#eta(MonoJet)", 50, -2.5, 2.5);
0200   monojet_PFJet_phi = bei.book1D("monojet_PFJet_phi", "#phi(MonoJet)", 50, -3.14, 3.14);
0201   monojet_PFMet = bei.book1D("monojet_PFMet", "Pt of PFMET (GeV)", 40, 0.0, 1000);
0202   monojet_PFMet_phi = bei.book1D("monojet_PFMet_phi", "#phi(PFMET #phi)", 50, -3.14, 3.14);
0203   monojet_PFJetPtOverPFMet = bei.book1D("monojet_PFJetPtOverPFMet", "Pt of MonoJet/MET (GeV)", 40, 0.0, 5.);
0204   monojet_deltaPhiPFJetPFMet = bei.book1D("monojet_deltaPhiPFJetPFMet", "#Delta#phi(MonoJet, PFMet)", 40, 0., 3.15);
0205   monojet_PFchef = bei.book1D("monojet_PFchef", "MonojetJet CHEF", 50, 0.0, 1.0);
0206   monojet_PFnhef = bei.book1D("monojet_PFnhef", "MonojetJet NHEF", 50, 0.0, 1.0);
0207   monojet_PFcemf = bei.book1D("monojet_PFcemf", "MonojetJet CEMF", 50, 0.0, 1.0);
0208   monojet_PFnemf = bei.book1D("monojet_PFnemf", "MonojetJet NEMF", 50, 0.0, 1.0);
0209   monojet_PFJetMulti = bei.book1D("monojet_PFJetMulti", "No. of PFJets", 10, 0., 10.);
0210   //--- MonoMuon
0211   bei.setCurrentFolder("Physics/Exotica/MonoMuon");
0212   monomuon_Muon_pt = bei.book1D("monomuon_Muon_pt", "Pt of Monomuon (GeV)", 50, 30.0, 2000);
0213   monomuon_Muon_eta = bei.book1D("monomuon_Muon_eta", "#eta(Monomuon)", 50, -2.5, 2.5);
0214   monomuon_Muon_phi = bei.book1D("monomuon_Muon_phi", "#phi(Monomuon)", 50, -3.14, 3.14);
0215   monomuon_Charge = bei.book1D("monomuon_Charge", "Charge of the MonoMuon", 10, -5., 5.);
0216   monomuon_PFMet = bei.book1D("monomuon_PFMet", "Pt of PFMET (GeV)", 40, 0.0, 2000);
0217   monomuon_PFMet_phi = bei.book1D("monomuon_PFMet_phi", "PFMET #phi", 50, -3.14, 3.14);
0218   monomuon_MuonPtOverPFMet = bei.book1D("monomuon_MuonPtOverPFMet", "Pt of Monomuon/PFMet", 40, 0.0, 5.);
0219   monomuon_deltaPhiMuonPFMet = bei.book1D("monomuon_deltaPhiMuonPFMet", "#Delta#phi(Monomuon, PFMet)", 40, 0., 3.15);
0220   monomuon_TransverseMass = bei.book1D("monomuon_TransverseMass", "Transverse Mass M_{T} GeV", 40, 200., 3000.);
0221   monomuon_MuonMulti = bei.book1D("monomuon_MuonMulti", "No. of Muons", 10, 0., 10.);
0222   //--- MonoElectron
0223   bei.setCurrentFolder("Physics/Exotica/MonoElectron");
0224   monoelectron_Electron_pt = bei.book1D("monoelectron_Electron_pt", "Pt of Monoelectron (GeV)", 50, 30.0, 4000);
0225   monoelectron_Electron_eta = bei.book1D("monoelectron_Electron_eta", "#eta(MonoElectron)", 50, -2.5, 2.5);
0226   monoelectron_Electron_phi = bei.book1D("monoelectron_Electron_phi", "#phi(MonoElectron)", 50, -3.14, 3.14);
0227   monoelectron_Charge = bei.book1D("monoelectron_Charge", "Charge of the MonoElectron", 10, -5., 5.);
0228   monoelectron_PFMet = bei.book1D("monoelectron_PFMet", "Pt of PFMET (GeV)", 40, 0.0, 4000);
0229   monoelectron_PFMet_phi = bei.book1D("monoelectron_PFMet_phi", "PFMET #phi", 50, -3.14, 3.14);
0230   monoelectron_ElectronPtOverPFMet =
0231       bei.book1D("monoelectron_ElectronPtOverPFMet", "Pt of Monoelectron/PFMet", 40, 0.0, 5.);
0232   monoelectron_deltaPhiElectronPFMet =
0233       bei.book1D("monoelectron_deltaPhiElectronPFMet", "#Delta#phi(MonoElectron, PFMet)", 40, 0., 3.15);
0234   monoelectron_TransverseMass = bei.book1D("monoelectron_TransverseMass", "Transverse Mass M_{T} GeV", 40, 200., 4000.);
0235   monoelectron_ElectronMulti = bei.book1D("monoelectron_ElectronMulti", "No. of Electrons", 10, 0., 10.);
0236 
0237   //--- DiPhotons
0238   bei.setCurrentFolder("Physics/Exotica/MonoPhotons");
0239   monophoton_Photon_energy = bei.book1D("monophoton_Photon_energy", "Energy of Leading Photon (GeV)", 50, 30.0, 1000);
0240   monophoton_Photon_et = bei.book1D("monophoton_Photon_et", "Et of Leading Photon (GeV)", 50, 30.0, 1000);
0241   monophoton_Photon_pt = bei.book1D("monophoton_Photon_pt", "Pt of Leading Photon (GeV)", 50, 30.0, 1000);
0242   monophoton_Photon_eta = bei.book1D("monophoton_Photon_eta", "#eta (Leading Photon)", 50, -2.5, 2.5);
0243   monophoton_Photon_etasc = bei.book1D("monophoton_Photon_etasc", "#eta sc(Leading Photon)", 50, -2.5, 2.5);
0244   monophoton_Photon_phi = bei.book1D("monophoton_Photon_phi", "#phi(Leading Photon)", 50, -3.14, 3.14);
0245   monophoton_Photon_hovere = bei.book1D("monophoton_Photon_hovere", "H/E (Leading Photon)", 50, 0., 0.50);
0246   monophoton_Photon_sigmaietaieta =
0247       bei.book1D("monophoton_Photon_sigmaietaieta", "#sigma_{i #eta i #eta} (Leading Photon)", 50, 0., 0.03);
0248   monophoton_Photon_trksumptsolidconedr03 =
0249       bei.book1D("monophoton_Photon_trksumptsolidconedr03", "TrkSumPtDr03 (Leading Photon)", 50, 0., 15.);
0250   monophoton_Photon_e1x5e5x5 = bei.book1D("monophoton_Photon_e1x5e5x5", "E_{1x5}/E_{5x5} (Leading Photon)", 50, 0., 1.);
0251   monophoton_Photon_e2x5e5x5 = bei.book1D("monophoton_Photon_e2x5e5x5", "E_{2x5}/E_{5x5} (Leading Photon)", 50, 0., 1.);
0252   monophoton_PFMet = bei.book1D("monophoton_PFMet", "Pt of PFMET (GeV)", 40, 0.0, 1000);
0253   monophoton_PFMet_phi = bei.book1D("monophoton_PFMet_phi", "PFMET #phi", 50, -3.14, 3.14);
0254   monophoton_PhotonPtOverPFMet = bei.book1D("monophoton_PhotonPtOverPFMet", "Pt of Monophoton/PFMet", 40, 0.0, 5.);
0255   monophoton_deltaPhiPhotonPFMet =
0256       bei.book1D("monophoton_deltaPhiPhotonPFMet", "#Delta#phi(SubLeading Photon, PFMet)", 40, 0., 3.15);
0257   monophoton_PhotonMulti = bei.book1D("monophoton_PhotonMulti", "No. of Photons", 10, 0., 10.);
0258 
0259   //--- Displaced Leptons (filled using only leptons from long-lived stop decay).
0260   bei.setCurrentFolder("Physics/Exotica/DisplacedFermions");
0261   dispElec_track_effi_lxy = bei.bookProfile("dispElec_track_effi_lxy",
0262                                             "Electron channel; transverse decay length (cm); track reco efficiency",
0263                                             10,
0264                                             0.,
0265                                             100.,
0266                                             -999.,
0267                                             999,
0268                                             "");
0269   dispElec_elec_effi_lxy = bei.bookProfile("dispElec_elec_effi_lxy",
0270                                            "Electron channel; transverse decay length (cm); electron reco efficiency",
0271                                            10,
0272                                            0.,
0273                                            100.,
0274                                            -999.,
0275                                            999,
0276                                            "");
0277   dispMuon_track_effi_lxy = bei.bookProfile("dispMuon_track_effi_lxy",
0278                                             "Muon channel; transverse decay length (cm); track reco efficiency",
0279                                             10,
0280                                             0.,
0281                                             100.,
0282                                             -999.,
0283                                             999,
0284                                             "");
0285   dispMuon_muon_effi_lxy = bei.bookProfile("dispMuon_muon_effi_lxy",
0286                                            "Muon channel; transverse decay length (cm); muon reco efficiency",
0287                                            10,
0288                                            0.,
0289                                            100.,
0290                                            -999.,
0291                                            999,
0292                                            "");
0293   dispMuon_muonDisp_effi_lxy =
0294       bei.bookProfile("dispMuon_muonDisp_effi_lxy",
0295                       "Muon channel; transverse decay length (cm); displacedMuon reco efficiency",
0296                       10,
0297                       0.,
0298                       100.,
0299                       -999.,
0300                       999,
0301                       "");
0302   dispMuon_muonDispSA_effi_lxy =
0303       bei.bookProfile("dispMuon_muonDispSA_effi_lxy",
0304                       "Muon channel; transverse decay length (cm); displacedSAMuon reco efficiency",
0305                       10,
0306                       0.,
0307                       400.,
0308                       -999.,
0309                       999,
0310                       "");
0311   //--- Displaced Jets (filled using only tracks or jets from long-lived stop decay).
0312   dispJet_track_effi_lxy = bei.bookProfile("dispJet_track_effi_lxy",
0313                                            "Jet channel; transverse decay length (cm); track reco efficiency",
0314                                            10,
0315                                            0.,
0316                                            100.,
0317                                            -999.,
0318                                            999,
0319                                            "");
0320 
0321   bei.cd();
0322 }
0323 
0324 //
0325 //  -- Analyze
0326 //
0327 void ExoticaDQM::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0328   // objects
0329 
0330   //Trigger
0331   bool ValidTriggers = iEvent.getByToken(TriggerToken_, TriggerResults_);
0332   if (!ValidTriggers)
0333     return;
0334 
0335   // Vertices
0336   bool ValidVertices = iEvent.getByToken(VertexToken_, VertexCollection_);
0337   if (!ValidVertices)
0338     return;
0339 
0340   // Electrons
0341   bool ValidGedGsfElectron = iEvent.getByToken(ElectronToken_, ElectronCollection_);
0342   if (!ValidGedGsfElectron)
0343     return;
0344 
0345   // Muons
0346   bool ValidPFMuon = iEvent.getByToken(MuonToken_, MuonCollection_);
0347   if (!ValidPFMuon)
0348     return;
0349 
0350   // Jets
0351   bool ValidPFJet = iEvent.getByToken(PFJetToken_, pfJetCollection_);
0352   pfjets = *pfJetCollection_;
0353   if (!ValidPFJet)
0354     return;
0355 
0356   // PFMETs
0357   bool ValidPFMET = iEvent.getByToken(PFMETToken_, pfMETCollection_);
0358   if (!ValidPFMET)
0359     return;
0360 
0361   // Photons
0362   bool ValidCaloPhoton = iEvent.getByToken(PhotonToken_, PhotonCollection_);
0363   if (!ValidCaloPhoton)
0364     return;
0365 
0366   // Tracks
0367   bool ValidTracks = iEvent.getByToken(TrackToken_, TrackCollection_);
0368   if (!ValidTracks)
0369     return;
0370 
0371   // Special collections for displaced particles
0372   iEvent.getByToken(MuonDispToken_, MuonDispCollection_);
0373   iEvent.getByToken(MuonDispSAToken_, MuonDispSACollection_);
0374 
0375   // MC truth
0376   bool ValidGenParticles = iEvent.getByToken(GenParticleToken_, GenCollection_);
0377 
0378   // JetCorrector
0379   bool ValidJetCorrector = iEvent.getByToken(JetCorrectorToken_, JetCorrector_);
0380 
0381   //Trigger
0382 
0383   int N_Triggers = TriggerResults_->size();
0384   int N_GoodTriggerPaths = HltPaths_.size();
0385   bool triggered_event = false;
0386   const edm::TriggerNames& trigName = iEvent.triggerNames(*TriggerResults_);
0387   for (int i_Trig = 0; i_Trig < N_Triggers; ++i_Trig) {
0388     if (TriggerResults_.product()->accept(i_Trig)) {
0389       for (int n = 0; n < N_GoodTriggerPaths; n++) {
0390         if (trigName.triggerName(i_Trig).find(HltPaths_[n]) != std::string::npos) {
0391           triggered_event = true;
0392         }
0393       }
0394     }
0395   }
0396   if (triggered_event == false)
0397     return;
0398 
0399   for (int i = 0; i < 2; i++) {
0400     //Jets
0401     PFJetPx[i] = 0.;
0402     PFJetPy[i] = 0.;
0403     PFJetPt[i] = 0.;
0404     PFJetEta[i] = 0.;
0405     PFJetPhi[i] = 0.;
0406     PFJetNHEF[i] = 0.;
0407     PFJetCHEF[i] = 0.;
0408     PFJetNEMF[i] = 0.;
0409     PFJetCEMF[i] = 0.;
0410     //Muons
0411     MuonPx[i] = 0.;
0412     MuonPy[i] = 0.;
0413     MuonPt[i] = 0.;
0414     MuonEta[i] = 0.;
0415     MuonPhi[i] = 0.;
0416     MuonCharge[i] = 0.;
0417     //Electrons
0418     ElectronPx[i] = 0.;
0419     ElectronPy[i] = 0.;
0420     ElectronPt[i] = 0.;
0421     ElectronEta[i] = 0.;
0422     ElectronPhi[i] = 0.;
0423     ElectronCharge[i] = 0.;
0424     //Photons
0425     PhotonEnergy[i] = 0.;
0426     PhotonPt[i] = 0.;
0427     PhotonEt[i] = 0.;
0428     PhotonEta[i] = 0.;
0429     PhotonEtaSc[i] = 0.;
0430     PhotonPhi[i] = 0.;
0431     PhotonHoverE[i] = 0.;
0432     PhotonSigmaIetaIeta[i] = 0.;
0433     PhotonTrkSumPtSolidConeDR03[i] = 0.;
0434     PhotonE1x5E5x5[i] = 0.;
0435     PhotonE2x5E5x5[i] = 0.;
0436   }
0437 
0438   //Getting information from the RecoObjects
0439   dijet_countPFJet_ = 0;
0440   monojet_countPFJet_ = 0;
0441 
0442   PFJetCollection::const_iterator pfjet_ = pfjets.begin();
0443   for (; pfjet_ != pfjets.end(); ++pfjet_) {
0444     double scale = 1.;
0445     if (ValidJetCorrector)
0446       scale = JetCorrector_->correction(*pfjet_);
0447     if (scale * pfjet_->pt() > PFJetPt[0]) {
0448       PFJetPt[1] = PFJetPt[0];
0449       PFJetPx[1] = PFJetPx[0];
0450       PFJetPy[1] = PFJetPy[0];
0451       PFJetEta[1] = PFJetEta[0];
0452       PFJetPhi[1] = PFJetPhi[0];
0453       PFJetRapidity[1] = PFJetRapidity[0];
0454       PFJetMass[1] = PFJetMass[0];
0455       PFJetNHEF[1] = PFJetNHEF[0];
0456       PFJetCHEF[1] = PFJetCHEF[0];
0457       PFJetNEMF[1] = PFJetNEMF[0];
0458       PFJetCEMF[1] = PFJetCEMF[0];
0459       //
0460       PFJetPt[0] = scale * pfjet_->pt();
0461       PFJetPx[0] = scale * pfjet_->px();
0462       PFJetPy[0] = scale * pfjet_->py();
0463       PFJetEta[0] = pfjet_->eta();
0464       PFJetPhi[0] = pfjet_->phi();
0465       PFJetRapidity[0] = pfjet_->rapidity();
0466       PFJetMass[0] = pfjet_->mass();
0467       PFJetNHEF[0] = pfjet_->neutralHadronEnergyFraction();
0468       PFJetCHEF[0] = pfjet_->chargedHadronEnergyFraction();
0469       PFJetNEMF[0] = pfjet_->neutralEmEnergyFraction();
0470       PFJetCEMF[0] = pfjet_->chargedEmEnergyFraction();
0471     } else if (scale * pfjet_->pt() < PFJetPt[0] && scale * pfjet_->pt() > PFJetPt[1]) {
0472       PFJetPt[1] = scale * pfjet_->pt();
0473       PFJetPx[1] = scale * pfjet_->px();
0474       PFJetPy[1] = scale * pfjet_->py();
0475       PFJetEta[1] = pfjet_->eta();
0476       PFJetPhi[1] = pfjet_->phi();
0477       PFJetRapidity[1] = pfjet_->rapidity();
0478       PFJetMass[1] = pfjet_->mass();
0479       PFJetNHEF[1] = pfjet_->neutralHadronEnergyFraction();
0480       PFJetCHEF[1] = pfjet_->chargedHadronEnergyFraction();
0481       PFJetNEMF[1] = pfjet_->neutralEmEnergyFraction();
0482       PFJetCEMF[1] = pfjet_->chargedEmEnergyFraction();
0483     } else {
0484     }
0485     if (scale * pfjet_->pt() > dijet_PFJet1_pt_cut_)
0486       dijet_countPFJet_++;
0487     if (scale * pfjet_->pt() > dijet_PFJet1_pt_cut_)
0488       monojet_countPFJet_++;
0489   }
0490 
0491   VertexCollection vertexCollection = *(VertexCollection_.product());
0492   reco::VertexCollection::const_iterator primaryVertex_ = vertexCollection.begin();
0493 
0494   dimuon_countMuon_ = 0;
0495   monomuon_countMuon_ = 0;
0496   reco::MuonCollection::const_iterator muon_ = MuonCollection_->begin();
0497   for (; muon_ != MuonCollection_->end(); muon_++) {
0498     // Muon High Pt ID
0499     bool HighPt = false;
0500     if (muon_->isGlobalMuon() && muon_->globalTrack()->hitPattern().numberOfValidMuonHits() > 0 &&
0501         muon_->numberOfMatchedStations() > 1 && muon_->innerTrack()->hitPattern().trackerLayersWithMeasurement() > 5 &&
0502         muon_->innerTrack()->hitPattern().numberOfValidPixelHits() > 0 &&
0503         muon_->muonBestTrack()->ptError() / muon_->muonBestTrack()->pt() < 0.3 &&
0504         fabs(muon_->muonBestTrack()->dxy(primaryVertex_->position())) < 0.2 &&
0505         fabs(muon_->bestTrack()->dz(primaryVertex_->position())) < 0.5 && fabs(muon_->eta()) < 2.1)
0506       HighPt = true;
0507 
0508     if (HighPt == true) {
0509       if (muon_->pt() > MuonPt[0]) {
0510         MuonPt[1] = MuonPt[0];
0511         MuonPx[1] = MuonPx[0];
0512         MuonPy[1] = MuonPy[0];
0513         MuonEta[1] = MuonEta[0];
0514         MuonPhi[1] = MuonPhi[0];
0515         MuonCharge[1] = MuonCharge[0];
0516         //
0517         MuonPt[0] = muon_->pt();
0518         MuonPx[0] = muon_->px();
0519         MuonPy[0] = muon_->py();
0520         MuonEta[0] = muon_->eta();
0521         MuonPhi[0] = muon_->phi();
0522         MuonCharge[0] = muon_->charge();
0523       } else if (muon_->pt() > MuonPt[1] && muon_->pt() < MuonPt[0]) {
0524         MuonPt[1] = muon_->pt();
0525         MuonPx[1] = muon_->px();
0526         MuonPy[1] = muon_->py();
0527         MuonEta[1] = muon_->eta();
0528         MuonPhi[1] = muon_->phi();
0529         MuonCharge[1] = muon_->charge();
0530       }
0531     }
0532     if (muon_->pt() > dimuon_Muon1_pt_cut_)
0533       dimuon_countMuon_++;
0534     if (muon_->pt() > monomuon_Muon_pt_cut_)
0535       monomuon_countMuon_++;
0536   }
0537 
0538   dielectron_countElectron_ = 0;
0539   monoelectron_countElectron_ = 0;
0540   reco::GsfElectronCollection::const_iterator electron_ = ElectronCollection_->begin();
0541   for (; electron_ != ElectronCollection_->end(); electron_++) {
0542     //HEEP Selection 4.1 (some cuts)
0543     if (electron_->e5x5() <= 0)
0544       continue;
0545     if (electron_->gsfTrack().isNull())
0546       continue;
0547     bool HEPP_ele = false;
0548     double sceta = electron_->caloPosition().eta();
0549     double dEtaIn = fabs(electron_->deltaEtaSuperClusterTrackAtVtx());
0550     double dPhiIn = fabs(electron_->deltaPhiSuperClusterTrackAtVtx());
0551     double HoverE = electron_->hadronicOverEm();
0552     int missingHits = electron_->gsfTrack()->hitPattern().numberOfLostTrackerHits(HitPattern::MISSING_INNER_HITS);
0553     double dxy = electron_->gsfTrack()->dxy(primaryVertex_->position());
0554     double tkIso = electron_->dr03TkSumPt();
0555     double e2x5Fraction = electron_->e2x5Max() / electron_->e5x5();
0556     double e1x5Fraction = electron_->e1x5() / electron_->e5x5();
0557     double scSigmaIetaIeta = electron_->scSigmaIEtaIEta();
0558     if (electron_->ecalDriven() && electron_->pt() > 35.) {
0559       if (fabs(sceta) < 1.442) {  // barrel
0560         if (fabs(dEtaIn) < 0.005 && fabs(dPhiIn) < 0.06 && HoverE < 0.05 && tkIso < 5. && missingHits <= 1 &&
0561             fabs(dxy) < 0.02 && (e2x5Fraction > 0.94 || e1x5Fraction > 0.83))
0562           HEPP_ele = true;
0563       } else if (fabs(sceta) > 1.56 && fabs(sceta) < 2.5) {  // endcap
0564         if (fabs(dEtaIn) < 0.007 && fabs(dPhiIn) < 0.06 && HoverE < 0.05 && tkIso < 5. && missingHits <= 1 &&
0565             fabs(dxy) < 0.02 && scSigmaIetaIeta < 0.03)
0566           HEPP_ele = true;
0567       }
0568     }
0569     //
0570     if (HEPP_ele == false)
0571       continue;
0572     if (electron_->pt() > ElectronPt[0]) {
0573       ElectronPt[1] = ElectronPt[0];
0574       ElectronPx[1] = ElectronPx[0];
0575       ElectronPy[1] = ElectronPy[0];
0576       ElectronEta[1] = ElectronEta[0];
0577       ElectronPhi[1] = ElectronPhi[0];
0578       ElectronCharge[1] = ElectronCharge[0];
0579       //
0580       ElectronPt[0] = electron_->pt();
0581       ElectronPx[0] = electron_->px();
0582       ElectronPy[0] = electron_->py();
0583       ElectronEta[0] = electron_->eta();
0584       ElectronPhi[0] = electron_->phi();
0585       ElectronCharge[0] = electron_->charge();
0586     } else if (electron_->pt() > ElectronPt[1] && electron_->pt() < ElectronPt[0]) {
0587       ElectronPt[1] = electron_->pt();
0588       ElectronPx[1] = electron_->px();
0589       ElectronPy[1] = electron_->py();
0590       ElectronEta[1] = electron_->eta();
0591       ElectronPhi[1] = electron_->phi();
0592       ElectronCharge[1] = electron_->charge();
0593     }
0594     if (electron_->pt() > dielectron_Electron1_pt_cut_)
0595       dielectron_countElectron_++;
0596     if (electron_->pt() > monoelectron_Electron_pt_cut_)
0597       monoelectron_countElectron_++;
0598   }
0599 
0600   diphoton_countPhoton_ = 0.;
0601   reco::PhotonCollection::const_iterator photon_ = PhotonCollection_->begin();
0602   for (; photon_ != PhotonCollection_->end(); ++photon_) {
0603     if (photon_->pt() > PhotonPt[0]) {
0604       PhotonEnergy[1] = PhotonEnergy[0];
0605       PhotonPt[1] = PhotonPt[0];
0606       PhotonEt[1] = PhotonEt[0];
0607       PhotonEta[1] = PhotonEta[0];
0608       PhotonEtaSc[1] = PhotonEtaSc[0];
0609       PhotonPhi[1] = PhotonPhi[0];
0610       PhotonHoverE[1] = PhotonHoverE[0];
0611       PhotonSigmaIetaIeta[1] = PhotonSigmaIetaIeta[0];
0612       PhotonTrkSumPtSolidConeDR03[1] = PhotonTrkSumPtSolidConeDR03[0];
0613       PhotonE1x5E5x5[1] = PhotonE1x5E5x5[0];
0614       PhotonE2x5E5x5[1] = PhotonE2x5E5x5[0];
0615       PhotonEnergy[0] = photon_->energy();
0616       PhotonPt[0] = photon_->pt();
0617       PhotonEt[0] = photon_->et();
0618       PhotonEta[0] = photon_->eta();
0619       PhotonEtaSc[0] = photon_->caloPosition().eta();
0620       PhotonPhi[0] = photon_->phi();
0621       PhotonHoverE[0] = photon_->hadronicOverEm();
0622       PhotonSigmaIetaIeta[0] = photon_->sigmaIetaIeta();
0623       PhotonTrkSumPtSolidConeDR03[0] = photon_->trkSumPtSolidConeDR03();
0624       PhotonE1x5E5x5[0] = photon_->e1x5() / photon_->e5x5();
0625       PhotonE2x5E5x5[0] = photon_->e2x5() / photon_->e5x5();
0626     } else if (photon_->pt() > PhotonPt[1] && photon_->pt() < PhotonPt[0]) {
0627       PhotonEnergy[1] = photon_->energy();
0628       PhotonPt[1] = photon_->pt();
0629       PhotonEt[1] = photon_->et();
0630       PhotonEta[1] = photon_->eta();
0631       PhotonEtaSc[1] = photon_->caloPosition().eta();
0632       PhotonPhi[1] = photon_->phi();
0633       PhotonHoverE[1] = photon_->hadronicOverEm();
0634       PhotonSigmaIetaIeta[1] = photon_->sigmaIetaIeta();
0635       PhotonTrkSumPtSolidConeDR03[1] = photon_->trkSumPtSolidConeDR03();
0636       PhotonE1x5E5x5[1] = photon_->e1x5() / photon_->e5x5();
0637       PhotonE2x5E5x5[1] = photon_->e2x5() / photon_->e5x5();
0638     }
0639     if (photon_->pt() > diphoton_Photon1_pt_cut_)
0640       diphoton_countPhoton_++;
0641   }
0642 
0643   //
0644   // Analyze
0645   //
0646 
0647   //Resonances
0648   analyzeDiJets(iEvent);
0649   analyzeDiMuons(iEvent);
0650   analyzeDiElectrons(iEvent);
0651   analyzeDiPhotons(iEvent);
0652 
0653   //MonoSearches
0654   analyzeMonoJets(iEvent);
0655   analyzeMonoMuons(iEvent);
0656   analyzeMonoElectrons(iEvent);
0657 
0658   //LongLived
0659   if (ValidGenParticles) {
0660     analyzeDisplacedLeptons(iEvent, iSetup);
0661     analyzeDisplacedJets(iEvent, iSetup);
0662   }
0663 }
0664 
0665 void ExoticaDQM::analyzeDisplacedLeptons(const Event& iEvent, const edm::EventSetup& iSetup) {
0666   //=== This is designed to run on MC events in which a pair of long-lived stop quarks each decay to a displaced lepton + displaced b jet.
0667 
0668   // Initialisation
0669 
0670   const unsigned int stop1 = 1000006;  // PDG identifier of top squark1
0671   const unsigned int stop2 = 2000006;  // PDG identifier of top squark2
0672   const float deltaRcut = 0.01;        // Cone size for matching reco to true leptons.
0673   const float invPtcut = 0.1;          // Cut in 1/Pt consistency for matching reco tracks to genParticles.
0674 
0675   //--- Measure the efficiency to reconstruct leptons from long-lived stop quark decay.
0676 
0677   for (const reco::GenParticle& gen : *GenCollection_) {
0678     unsigned int idPdg = abs(gen.pdgId());
0679     // Find electrons/muons from long-lived stop decay.
0680     if (idPdg == stop1 || idPdg == stop2) {
0681       unsigned int nDau = gen.numberOfDaughters();
0682       for (unsigned int i = 0; i < nDau; i++) {
0683         const reco::GenParticle* dau = (const reco::GenParticle*)gen.daughter(i);
0684         // Only measure efficiency using leptons passing pt & eta cuts. (The pt cut is almost irrelevant, since leptons from stop decay are hard).
0685         if (fabs(dau->eta()) < dispFermion_eta_cut_ && dau->pt() > dispFermion_pt_cut_) {
0686           unsigned int pdgIdDau = abs(dau->pdgId());
0687 
0688           if (pdgIdDau == 11 || pdgIdDau == 13) {  // electron or muon from stop decay
0689 
0690             // Get transverse decay length of stop quark.
0691             float lxy = dau->vertex().rho();
0692 
0693             // Get momentum vector of daughter genParticle trajectory extrapolated to beam-line.
0694             GlobalVector genP = this->getGenParticleTrajectoryAtBeamline(iSetup, dau);
0695 
0696             if (pdgIdDau == 11) {  // electron from stop decay
0697 
0698               // Find matching reco track if any.
0699               bool recoedTrk = false;
0700               for (const reco::Track& trk : *TrackCollection_) {
0701                 if (reco::deltaR(genP, trk) < deltaRcut && fabs(1 / dau->pt() - 1 / trk.pt()) < invPtcut) {
0702                   //cout<<"MATCH ELEC TRK "<<dau->pt()<<" "<<trk.pt()<<" "<<reco::deltaR(genP, trk)<<endl;
0703                   recoedTrk = true;
0704                 }
0705               }
0706               dispElec_track_effi_lxy->Fill(lxy, recoedTrk);
0707 
0708               // Find matching reco electron if any.
0709               bool recoedE = false;
0710               for (const reco::GsfElectron& eReco : *ElectronCollection_) {
0711                 if (reco::deltaR(genP, eReco) < deltaRcut && fabs(1 / dau->pt() - 1 / eReco.pt()) < invPtcut)
0712                   recoedE = true;
0713               }
0714               dispElec_elec_effi_lxy->Fill(lxy, recoedE);
0715 
0716             } else if (pdgIdDau == 13) {  // muon from stop decay
0717 
0718               // Find matching reco track if any.
0719               bool recoedTrk = false;
0720               for (const reco::Track& trk : *TrackCollection_) {
0721                 if (reco::deltaR(genP, trk) < deltaRcut && fabs(1 / dau->pt() - 1 / trk.pt()) < invPtcut) {
0722                   //cout<<"MATCH MUON TRK "<<dau->pt()<<" "<<trk.pt()<<" "<<reco::deltaR(genP, trk)<<endl;
0723                   recoedTrk = true;
0724                 }
0725               }
0726               dispMuon_track_effi_lxy->Fill(lxy, recoedTrk);
0727 
0728               // Find matching reco muon, if any, in normal global muon collection.
0729               bool recoedMu = false;
0730               for (const reco::Muon& muReco : *MuonCollection_) {
0731                 if (reco::deltaR(genP, muReco) < deltaRcut && fabs(1 / dau->pt() - 1 / muReco.pt()) < invPtcut)
0732                   recoedMu = true;
0733               }
0734               dispMuon_muon_effi_lxy->Fill(lxy, recoedMu);
0735 
0736               // Find matching reco muon, if any, in displaced global muon collection.
0737               bool recoedMuDisp = false;
0738               for (const reco::Track& muDispReco : *MuonDispCollection_) {
0739                 if (reco::deltaR(genP, muDispReco) < deltaRcut && fabs(1 / dau->pt() - 1 / muDispReco.pt()) < invPtcut)
0740                   recoedMuDisp = true;
0741               }
0742               dispMuon_muonDisp_effi_lxy->Fill(lxy, recoedMuDisp);
0743 
0744               // Find matching reco muon, if any, in displaced SA muon collection.
0745               bool recoedMuDispSA = false;
0746               for (const reco::Track& muDispSAReco : *MuonDispSACollection_) {
0747                 if (reco::deltaR(genP, muDispSAReco) < deltaRcut &&
0748                     fabs(1 / dau->pt() - 1 / muDispSAReco.pt()) < invPtcut)
0749                   recoedMuDispSA = true;
0750               }
0751               dispMuon_muonDispSA_effi_lxy->Fill(lxy, recoedMuDispSA);
0752             }
0753           }
0754         }
0755       }
0756     }
0757   }
0758 }
0759 void ExoticaDQM::analyzeDisplacedJets(const Event& iEvent, const edm::EventSetup& iSetup) {
0760   //=== This is designed to run on MC events in which a pair of long-lived stop quarks each decay to a displaced lepton + displaced b jet.
0761 
0762   // Initialisation
0763 
0764   // Define function to identify R-hadrons containing stop quarks from PDG particle code.
0765   // N.B. Jets originate not just from stop quark, but also from its partner SM quark inside the R hadron.
0766   auto isRhadron = [](unsigned int pdgId) { return (pdgId / 100) == 10006 || (pdgId / 1000) == 1006; };
0767 
0768   const float deltaRcut = 0.01;  // Cone size for matching reco tracks to genParticles.
0769   const float invPtcut = 0.1;    // Cut in 1/Pt consistency for matching reco tracks to genParticles.
0770 
0771   //--- Measure the efficiency to reconstruct tracks in jet(s) from long-lived stop quark decay.
0772 
0773   for (const reco::GenParticle& gen : *GenCollection_) {
0774     unsigned int idPdg = abs(gen.pdgId());
0775     // Only measure efficiency using charged e, mu pi, K, p
0776     if (idPdg == 11 || idPdg == 13 || idPdg == 211 || idPdg == 321 || idPdg == 2212) {
0777       // Only measure efficiency using leptons passing pt & eta cuts. (The pt cut is almost irrelevant, since leptons from stop decay are hard).
0778       if (fabs(gen.eta()) < dispFermion_eta_cut_ && gen.pt() > dispFermion_pt_cut_) {
0779         // Check if this particle came (maybe indirectly) from an R hadron decay.
0780         const reco::GenParticle* genMoth = &gen;
0781         const reco::GenParticle* genRhadron = nullptr;
0782         bool foundParton = false;
0783         while (genMoth->numberOfMothers() > 0) {
0784           genMoth = (const reco::GenParticle*)genMoth->mother(0);
0785           unsigned int idPdgMoth = abs(genMoth->pdgId());
0786           // Check that the R-hadron decayed via a quark/gluon before yielding genParticle "gen".
0787           // This ensures that gen is from the jet, and not a lepton produced directly from the stop quark decay.
0788           if ((idPdgMoth >= 1 && idPdgMoth <= 6) || idPdgMoth == 21)
0789             foundParton = true;
0790           // Note if ancestor was R hadron
0791           if (isRhadron(idPdgMoth)) {
0792             genRhadron = genMoth;
0793             break;
0794           }
0795         }
0796 
0797         if (foundParton && genRhadron != nullptr) {  // This GenParticle came (maybe indirectly) from an R hadron decay.
0798 
0799           // Get transverse decay length of R hadron.
0800           float lxy = genRhadron->daughter(0)->vertex().rho();
0801 
0802           // Get momentum vector of genParticle trajectory extrapolated to beam-line.
0803           GlobalVector genP = this->getGenParticleTrajectoryAtBeamline(iSetup, &gen);
0804 
0805           // Find matching reco track if any.
0806           bool recoedTrk = false;
0807           for (const reco::Track& trk : *TrackCollection_) {
0808             if (reco::deltaR(genP, trk) < deltaRcut && fabs(1 / gen.pt() - 1 / trk.pt()) < invPtcut) {
0809               //cout<<"MATCH TRK "<<gen.pt()<<" "<<trk.pt()<<" "<<reco::deltaR(gen, trk)<<endl;
0810               recoedTrk = true;
0811             }
0812           }
0813           dispJet_track_effi_lxy->Fill(lxy, recoedTrk);
0814         }
0815       }
0816     }
0817   }
0818 }
0819 GlobalVector ExoticaDQM::getGenParticleTrajectoryAtBeamline(const edm::EventSetup& iSetup,
0820                                                             const reco::GenParticle* gen) {
0821   //=== Estimate the momentum vector that a GenParticle would have at its trajectory's point of closest
0822   //=== approach to the beam-line.
0823 
0824   // Get the magnetic field
0825   const MagneticField* theMagField = &iSetup.getData(magFieldToken_);
0826 
0827   // Make FreeTrajectoryState of this gen particle
0828   FreeTrajectoryState fts(GlobalPoint(gen->vx(), gen->vy(), gen->vz()),
0829                           GlobalVector(gen->px(), gen->py(), gen->pz()),
0830                           gen->charge(),
0831                           theMagField);
0832 
0833   // Get trajectory closest to beam line
0834   TSCBLBuilderNoMaterial tscblBuilder;
0835   const BeamSpot beamspot;  // Simple beam-spot at (0,0,0). Good enough.
0836   TrajectoryStateClosestToBeamLine tsAtClosestApproach = tscblBuilder(fts, beamspot);
0837 
0838   GlobalVector p = tsAtClosestApproach.trackStateAtPCA().momentum();
0839 
0840   return p;
0841 }
0842 
0843 void ExoticaDQM::analyzeDiJets(const Event& iEvent) {
0844   for (unsigned int icoll = 0; icoll < DiJetPFJetCollection_.size(); ++icoll) {
0845     dijet_countPFJet_ = 0;
0846     bool ValidDiJetPFJets = iEvent.getByToken(DiJetPFJetToken_[icoll], DiJetpfJetCollection_);
0847     if (!ValidDiJetPFJets)
0848       continue;
0849     DiJetpfjets = *DiJetpfJetCollection_;
0850     for (int i = 0; i < 2; i++) {
0851       PFJetPx[i] = 0.;
0852       PFJetPy[i] = 0.;
0853       PFJetPt[i] = 0.;
0854       PFJetEta[i] = 0.;
0855       PFJetPhi[i] = 0.;
0856       PFJetNHEF[i] = 0.;
0857       PFJetCHEF[i] = 0.;
0858       PFJetNEMF[i] = 0.;
0859       PFJetCEMF[i] = 0.;
0860     }
0861     PFJetCollection::const_iterator DiJetpfjet_ = DiJetpfjets.begin();
0862     for (; DiJetpfjet_ != DiJetpfjets.end(); ++DiJetpfjet_) {
0863       double scale = 1.;
0864       if (scale * DiJetpfjet_->pt() > PFJetPt[0]) {
0865         PFJetPt[1] = PFJetPt[0];
0866         PFJetPx[1] = PFJetPx[0];
0867         PFJetPy[1] = PFJetPy[0];
0868         PFJetEta[1] = PFJetEta[0];
0869         PFJetPhi[1] = PFJetPhi[0];
0870         PFJetRapidity[1] = DiJetpfjet_->rapidity();
0871         PFJetMass[1] = DiJetpfjet_->mass();
0872         PFJetNHEF[1] = PFJetNHEF[0];
0873         PFJetCHEF[1] = PFJetCHEF[0];
0874         PFJetNEMF[1] = PFJetNEMF[0];
0875         PFJetCEMF[1] = PFJetCEMF[0];
0876         PFJetPt[0] = scale * DiJetpfjet_->pt();
0877         PFJetPx[0] = scale * DiJetpfjet_->px();
0878         PFJetPy[0] = scale * DiJetpfjet_->py();
0879         PFJetEta[0] = DiJetpfjet_->eta();
0880         PFJetPhi[0] = DiJetpfjet_->phi();
0881         PFJetRapidity[0] = DiJetpfjet_->rapidity();
0882         PFJetMass[0] = DiJetpfjet_->mass();
0883         PFJetNHEF[0] = DiJetpfjet_->neutralHadronEnergyFraction();
0884         PFJetCHEF[0] = DiJetpfjet_->chargedHadronEnergyFraction();
0885         PFJetNEMF[0] = DiJetpfjet_->neutralEmEnergyFraction();
0886         PFJetCEMF[0] = DiJetpfjet_->chargedEmEnergyFraction();
0887       } else if (scale * DiJetpfjet_->pt() < PFJetPt[0] && scale * DiJetpfjet_->pt() > PFJetPt[1]) {
0888         PFJetPt[1] = scale * DiJetpfjet_->pt();
0889         PFJetPx[1] = scale * DiJetpfjet_->px();
0890         PFJetPy[1] = scale * DiJetpfjet_->py();
0891         PFJetEta[1] = DiJetpfjet_->eta();
0892         PFJetPhi[1] = DiJetpfjet_->phi();
0893         PFJetRapidity[1] = DiJetpfjet_->rapidity();
0894         PFJetMass[1] = DiJetpfjet_->mass();
0895         PFJetNHEF[1] = DiJetpfjet_->neutralHadronEnergyFraction();
0896         PFJetCHEF[1] = DiJetpfjet_->chargedHadronEnergyFraction();
0897         PFJetNEMF[1] = DiJetpfjet_->neutralEmEnergyFraction();
0898         PFJetCEMF[1] = DiJetpfjet_->chargedEmEnergyFraction();
0899       } else {
0900       }
0901       if (scale * DiJetpfjet_->pt() > dijet_PFJet1_pt_cut_)
0902         dijet_countPFJet_++;
0903     }
0904     if (PFJetPt[0] > dijet_PFJet1_pt_cut_ && PFJetPt[1] > dijet_PFJet2_pt_cut_) {
0905       dijet_PFJet_pt[icoll]->Fill(PFJetPt[0]);
0906       dijet_PFJet_eta[icoll]->Fill(PFJetEta[0]);
0907       dijet_PFJet_phi[icoll]->Fill(PFJetPhi[0]);
0908       dijet_PFJet_rapidity[icoll]->Fill(PFJetRapidity[0]);
0909       dijet_PFJet_mass[icoll]->Fill(PFJetMass[0]);
0910       dijet_PFJet_pt[icoll]->Fill(PFJetPt[1]);
0911       dijet_PFJet_eta[icoll]->Fill(PFJetEta[1]);
0912       dijet_PFJet_phi[icoll]->Fill(PFJetPhi[1]);
0913       dijet_PFJet_rapidity[icoll]->Fill(PFJetRapidity[1]);
0914       dijet_PFJet_mass[icoll]->Fill(PFJetMass[1]);
0915       dijet_deltaPhiPFJet1PFJet2[icoll]->Fill(deltaPhi(PFJetPhi[0], PFJetPhi[1]));
0916       dijet_deltaEtaPFJet1PFJet2[icoll]->Fill(PFJetEta[0] - PFJetEta[1]);
0917       dijet_deltaRPFJet1PFJet2[icoll]->Fill(deltaR(PFJetEta[0], PFJetPhi[0], PFJetEta[1], PFJetPhi[1]));
0918       dijet_invMassPFJet1PFJet2[icoll]->Fill(
0919           sqrt(2 * PFJetPt[0] * PFJetPt[1] * (cosh(PFJetEta[0] - PFJetEta[1]) - cos(PFJetPhi[0] - PFJetPhi[1]))));
0920       dijet_PFchef[icoll]->Fill(PFJetCHEF[0]);
0921       dijet_PFnhef[icoll]->Fill(PFJetNHEF[0]);
0922       dijet_PFcemf[icoll]->Fill(PFJetCEMF[0]);
0923       dijet_PFnemf[icoll]->Fill(PFJetNEMF[0]);
0924       dijet_PFJetMulti[icoll]->Fill(dijet_countPFJet_);
0925     }
0926   }
0927 }
0928 
0929 void ExoticaDQM::analyzeDiMuons(const Event& iEvent) {
0930   if (MuonPt[0] > dimuon_Muon1_pt_cut_ && MuonPt[1] > dimuon_Muon2_pt_cut_ && MuonCharge[0] * MuonCharge[1] == -1) {
0931     dimuon_Muon_pt->Fill(MuonPt[0]);
0932     dimuon_Muon_eta->Fill(MuonEta[0]);
0933     dimuon_Muon_phi->Fill(MuonPhi[0]);
0934     dimuon_Muon_pt->Fill(MuonPt[1]);
0935     dimuon_Muon_eta->Fill(MuonEta[1]);
0936     dimuon_Muon_phi->Fill(MuonPhi[1]);
0937     dimuon_Charge->Fill(MuonCharge[0]);
0938     dimuon_Charge->Fill(MuonCharge[1]);
0939     dimuon_deltaPhiMuon1Muon2->Fill(deltaPhi(MuonPhi[0], MuonPhi[1]));
0940     dimuon_deltaEtaMuon1Muon2->Fill(MuonEta[0] - MuonEta[1]);
0941     dimuon_deltaRMuon1Muon2->Fill(deltaR(MuonEta[0], MuonPhi[0], MuonEta[1], MuonPhi[1]));
0942     dimuon_invMassMuon1Muon2->Fill(
0943         sqrt(2 * MuonPt[0] * MuonPt[1] * (cosh(MuonEta[0] - MuonEta[1]) - cos(MuonPhi[0] - MuonPhi[1]))));
0944     dimuon_MuonMulti->Fill(dimuon_countMuon_);
0945   }
0946 }
0947 
0948 void ExoticaDQM::analyzeDiElectrons(const Event& iEvent) {
0949   if (ElectronPt[0] > dielectron_Electron1_pt_cut_ && ElectronPt[1] > dielectron_Electron2_pt_cut_ &&
0950       ElectronCharge[0] * ElectronCharge[1] == -1.) {
0951     dielectron_Electron_pt->Fill(ElectronPt[0]);
0952     dielectron_Electron_eta->Fill(ElectronEta[0]);
0953     dielectron_Electron_phi->Fill(ElectronPhi[0]);
0954     dielectron_Electron_pt->Fill(ElectronPt[1]);
0955     dielectron_Electron_eta->Fill(ElectronEta[1]);
0956     dielectron_Electron_phi->Fill(ElectronPhi[1]);
0957     dielectron_Charge->Fill(ElectronCharge[0]);
0958     dielectron_Charge->Fill(ElectronCharge[1]);
0959     dielectron_deltaPhiElectron1Electron2->Fill(deltaPhi(ElectronPhi[0], ElectronPhi[1]));
0960     dielectron_deltaEtaElectron1Electron2->Fill(ElectronEta[0] - ElectronEta[1]);
0961     dielectron_deltaRElectron1Electron2->Fill(deltaR(ElectronEta[0], ElectronPhi[0], ElectronEta[1], ElectronPhi[1]));
0962     dielectron_invMassElectron1Electron2->Fill(
0963         sqrt(2 * ElectronPt[0] * ElectronPt[1] *
0964              (cosh(ElectronEta[0] - ElectronEta[1]) - cos(ElectronPhi[0] - ElectronPhi[1]))));
0965     dielectron_ElectronMulti->Fill(dielectron_countElectron_);
0966   }
0967 }
0968 
0969 void ExoticaDQM::analyzeDiPhotons(const Event& iEvent) {
0970   if (PhotonPt[0] > diphoton_Photon1_pt_cut_ && PhotonPt[1] > diphoton_Photon2_pt_cut_) {
0971     diphoton_Photon_energy->Fill(PhotonEnergy[0]);
0972     diphoton_Photon_pt->Fill(PhotonPt[0]);
0973     diphoton_Photon_et->Fill(PhotonEt[0]);
0974     diphoton_Photon_eta->Fill(PhotonEta[0]);
0975     diphoton_Photon_etasc->Fill(PhotonEtaSc[0]);
0976     diphoton_Photon_phi->Fill(PhotonPhi[0]);
0977     if (fabs(PhotonEtaSc[0]) < 1.442) {
0978       diphoton_Photon_hovere_eb->Fill(PhotonHoverE[0]);
0979       diphoton_Photon_sigmaietaieta_eb->Fill(PhotonSigmaIetaIeta[0]);
0980       diphoton_Photon_trksumptsolidconedr03_eb->Fill(PhotonTrkSumPtSolidConeDR03[0]);
0981       diphoton_Photon_e1x5e5x5_eb->Fill(PhotonE1x5E5x5[0]);
0982       diphoton_Photon_e2x5e5x5_eb->Fill(PhotonE2x5E5x5[0]);
0983     }
0984     if (fabs(PhotonEtaSc[0]) > 1.566 && fabs(PhotonEtaSc[0]) < 2.5) {
0985       diphoton_Photon_hovere_ee->Fill(PhotonHoverE[0]);
0986       diphoton_Photon_sigmaietaieta_ee->Fill(PhotonSigmaIetaIeta[0]);
0987       diphoton_Photon_trksumptsolidconedr03_ee->Fill(PhotonTrkSumPtSolidConeDR03[0]);
0988       diphoton_Photon_e1x5e5x5_ee->Fill(PhotonE1x5E5x5[0]);
0989       diphoton_Photon_e2x5e5x5_ee->Fill(PhotonE2x5E5x5[0]);
0990     }
0991     diphoton_Photon_energy->Fill(PhotonEnergy[1]);
0992     diphoton_Photon_pt->Fill(PhotonPt[1]);
0993     diphoton_Photon_et->Fill(PhotonEt[1]);
0994     diphoton_Photon_eta->Fill(PhotonEta[1]);
0995     diphoton_Photon_etasc->Fill(PhotonEtaSc[1]);
0996     diphoton_Photon_phi->Fill(PhotonPhi[1]);
0997     if (fabs(PhotonEtaSc[1]) < 1.4442) {
0998       diphoton_Photon_hovere_eb->Fill(PhotonHoverE[1]);
0999       diphoton_Photon_sigmaietaieta_eb->Fill(PhotonSigmaIetaIeta[1]);
1000       diphoton_Photon_trksumptsolidconedr03_eb->Fill(PhotonTrkSumPtSolidConeDR03[1]);
1001       diphoton_Photon_e1x5e5x5_eb->Fill(PhotonE1x5E5x5[1]);
1002       diphoton_Photon_e2x5e5x5_eb->Fill(PhotonE2x5E5x5[1]);
1003     }
1004     if (fabs(PhotonEtaSc[1]) > 1.566 && fabs(PhotonEtaSc[1]) < 2.5) {
1005       diphoton_Photon_hovere_ee->Fill(PhotonHoverE[1]);
1006       diphoton_Photon_sigmaietaieta_ee->Fill(PhotonSigmaIetaIeta[1]);
1007       diphoton_Photon_trksumptsolidconedr03_ee->Fill(PhotonTrkSumPtSolidConeDR03[1]);
1008       diphoton_Photon_e1x5e5x5_ee->Fill(PhotonE1x5E5x5[1]);
1009       diphoton_Photon_e2x5e5x5_ee->Fill(PhotonE2x5E5x5[1]);
1010     }
1011     diphoton_deltaPhiPhoton1Photon2->Fill(deltaPhi(PhotonPhi[0], PhotonPhi[1]));
1012     diphoton_deltaEtaPhoton1Photon2->Fill(PhotonEta[0] - PhotonEta[1]);
1013     diphoton_deltaRPhoton1Photon2->Fill(deltaR(PhotonEta[0], PhotonPhi[0], PhotonEta[1], PhotonPhi[1]));
1014     diphoton_invMassPhoton1Photon2->Fill(
1015         sqrt(2 * PhotonPt[0] * PhotonPt[1] * (cosh(PhotonEta[0] - PhotonEta[1]) - cos(PhotonPhi[0] - PhotonPhi[1]))));
1016     diphoton_PhotonMulti->Fill(diphoton_countPhoton_);
1017   }
1018 }
1019 
1020 void ExoticaDQM::analyzeMonoJets(const Event& iEvent) {
1021   const PFMETCollection* pfmetcol = pfMETCollection_.product();
1022   const PFMET pfmet = pfmetcol->front();
1023   if (PFJetPt[0] > monojet_PFJet_pt_cut_ && pfmet.et() > monojet_PFJet_met_cut_) {
1024     monojet_PFJet_pt->Fill(PFJetPt[0]);
1025     monojet_PFJet_eta->Fill(PFJetEta[0]);
1026     monojet_PFJet_phi->Fill(PFJetPhi[0]);
1027     monojet_PFMet->Fill(pfmet.et());
1028     monojet_PFMet_phi->Fill(pfmet.phi());
1029     monojet_PFJetPtOverPFMet->Fill(PFJetPt[0] / pfmet.et());
1030     monojet_deltaPhiPFJetPFMet->Fill(deltaPhi(PFJetPhi[0], pfmet.phi()));
1031     monojet_PFchef->Fill(PFJetCHEF[0]);
1032     monojet_PFnhef->Fill(PFJetNHEF[0]);
1033     monojet_PFcemf->Fill(PFJetCEMF[0]);
1034     monojet_PFnemf->Fill(PFJetNEMF[0]);
1035     monojet_PFJetMulti->Fill(monojet_countPFJet_);
1036   }
1037 }
1038 
1039 void ExoticaDQM::analyzeMonoMuons(const Event& iEvent) {
1040   const PFMETCollection* pfmetcol = pfMETCollection_.product();
1041   const PFMET pfmet = pfmetcol->front();
1042   if (MuonPt[0] > monomuon_Muon_pt_cut_ && pfmet.et() > monomuon_Muon_met_cut_) {
1043     monomuon_Muon_pt->Fill(MuonPt[0]);
1044     monomuon_Muon_eta->Fill(MuonEta[0]);
1045     monomuon_Muon_phi->Fill(MuonPhi[0]);
1046     monomuon_Charge->Fill(MuonCharge[0]);
1047     monomuon_PFMet->Fill(pfmet.et());
1048     monomuon_PFMet_phi->Fill(pfmet.phi());
1049     monomuon_MuonPtOverPFMet->Fill(MuonPt[0] / pfmet.et());
1050     monomuon_deltaPhiMuonPFMet->Fill(deltaPhi(MuonPhi[0], pfmet.phi()));
1051     monomuon_TransverseMass->Fill(sqrt(2 * MuonPt[0] * pfmet.et() * (1 - cos(deltaPhi(MuonPhi[0], pfmet.phi())))));
1052     monomuon_MuonMulti->Fill(monomuon_countMuon_);
1053   }
1054 }
1055 
1056 void ExoticaDQM::analyzeMonoElectrons(const Event& iEvent) {
1057   const PFMETCollection* pfmetcol = pfMETCollection_.product();
1058   const PFMET pfmet = pfmetcol->front();
1059   if (ElectronPt[0] > monoelectron_Electron_pt_cut_ && pfmet.et() > monoelectron_Electron_met_cut_) {
1060     monoelectron_Electron_pt->Fill(ElectronPt[0]);
1061     monoelectron_Electron_eta->Fill(ElectronEta[0]);
1062     monoelectron_Electron_phi->Fill(ElectronPhi[0]);
1063     monoelectron_Charge->Fill(ElectronCharge[0]);
1064     monoelectron_PFMet->Fill(pfmet.et());
1065     monoelectron_PFMet_phi->Fill(pfmet.phi());
1066     monoelectron_ElectronPtOverPFMet->Fill(ElectronPt[0] / pfmet.et());
1067     monoelectron_deltaPhiElectronPFMet->Fill(deltaPhi(ElectronPhi[0], pfmet.phi()));
1068     monoelectron_TransverseMass->Fill(
1069         sqrt(2 * ElectronPt[0] * pfmet.et() * (1 - cos(deltaPhi(ElectronPhi[0], pfmet.phi())))));
1070     monoelectron_ElectronMulti->Fill(monoelectron_countElectron_);
1071   }
1072 }
1073 
1074 void ExoticaDQM::analyzeMonoPhotons(const Event& iEvent) {
1075   const PFMETCollection* pfmetcol = pfMETCollection_.product();
1076   const PFMET pfmet = pfmetcol->front();
1077   if (PhotonPt[0] > monophoton_Photon_pt_cut_ && pfmet.et() > monophoton_Photon_met_cut_) {
1078     monophoton_Photon_energy->Fill(PhotonEnergy[0]);
1079     monophoton_Photon_pt->Fill(PhotonPt[0]);
1080     monophoton_Photon_et->Fill(PhotonEt[0]);
1081     monophoton_Photon_eta->Fill(PhotonEta[0]);
1082     monophoton_Photon_etasc->Fill(PhotonEtaSc[0]);
1083     monophoton_Photon_phi->Fill(PhotonPhi[0]);
1084     monophoton_Photon_hovere->Fill(PhotonHoverE[0]);
1085     monophoton_Photon_sigmaietaieta->Fill(PhotonSigmaIetaIeta[0]);
1086     monophoton_Photon_trksumptsolidconedr03->Fill(PhotonTrkSumPtSolidConeDR03[0]);
1087     monophoton_Photon_e1x5e5x5->Fill(PhotonE1x5E5x5[0]);
1088     monophoton_Photon_e2x5e5x5->Fill(PhotonE2x5E5x5[0]);
1089     monophoton_PFMet->Fill(pfmet.et());
1090     monophoton_PFMet_phi->Fill(pfmet.phi());
1091     monophoton_PhotonPtOverPFMet->Fill(PhotonPt[0] / pfmet.et());
1092     monophoton_deltaPhiPhotonPFMet->Fill(deltaPhi(PhotonPhi[0], pfmet.phi()));
1093     monophoton_PhotonMulti->Fill(monophoton_countPhoton_);
1094   }
1095 }