00001 #include <cmath>
00002 #include <algorithm>
00003 #include <utility>
00004 #include <boost/foreach.hpp>
00005
00006 #include <TTree.h>
00007
00008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00009 #include "FWCore/Framework/interface/ESHandle.h"
00010 #include "FWCore/Framework/interface/EventSetup.h"
00011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00012
00013 #include "DataFormats/Common/interface/View.h"
00014 #include "DataFormats/JetReco/interface/Jet.h"
00015 #include "DataFormats/BTauReco/interface/JetTag.h"
00016
00017 #include "HLTrigger/HLTanalyzers/interface/HLTBJet.h"
00018
00019 static const size_t kMaxBJets = 10;
00020
00021 #include "HLTMessages.h"
00022
00023 HLTBJet::HLTBJet()
00024 {
00025
00026 NohBJetL2 = 0;
00027 ohBJetL2Energy = new float[kMaxBJets];
00028 ohBJetL2Et = new float[kMaxBJets];
00029 ohBJetL2Pt = new float[kMaxBJets];
00030 ohBJetL2Eta = new float[kMaxBJets];
00031 ohBJetL2Phi = new float[kMaxBJets];
00032
00033
00034 NohBJetL2Corrected = 0;
00035 ohBJetL2CorrectedEnergy = new float[kMaxBJets];
00036 ohBJetL2CorrectedEt = new float[kMaxBJets];
00037 ohBJetL2CorrectedPt = new float[kMaxBJets];
00038 ohBJetL2CorrectedEta = new float[kMaxBJets];
00039 ohBJetL2CorrectedPhi = new float[kMaxBJets];
00040
00041
00042 ohBJetIPL25Tag = new float[kMaxBJets];
00043 ohBJetIPL3Tag = new float[kMaxBJets];
00044
00045
00046 ohBJetIPLooseL25Tag = new float[kMaxBJets];
00047 ohBJetIPLooseL3Tag = new float[kMaxBJets];
00048
00049
00050 ohBJetMuL25Tag = new int[kMaxBJets];
00051 ohBJetMuL3Tag = new float[kMaxBJets];
00052
00053
00054 ohBJetPerfL25Tag = new int[kMaxBJets];
00055 ohBJetPerfL3Tag = new int[kMaxBJets];
00056 }
00057
00058 void HLTBJet::clear()
00059 {
00060 NohBJetL2 = 0;
00061 NohBJetL2Corrected = 0;
00062 std::memset(ohBJetL2Energy, '\0', kMaxBJets * sizeof(float));
00063 std::memset(ohBJetL2Et, '\0', kMaxBJets * sizeof(float));
00064 std::memset(ohBJetL2Et, '\0', kMaxBJets * sizeof(float));
00065 std::memset(ohBJetL2Pt, '\0', kMaxBJets * sizeof(float));
00066 std::memset(ohBJetL2Eta, '\0', kMaxBJets * sizeof(float));
00067 std::memset(ohBJetL2Phi, '\0', kMaxBJets * sizeof(float));
00068 std::memset(ohBJetL2CorrectedEnergy, '\0', kMaxBJets * sizeof(float));
00069 std::memset(ohBJetL2CorrectedEt, '\0', kMaxBJets * sizeof(float));
00070 std::memset(ohBJetL2CorrectedPt, '\0', kMaxBJets * sizeof(float));
00071 std::memset(ohBJetL2CorrectedEta, '\0', kMaxBJets * sizeof(float));
00072 std::memset(ohBJetL2CorrectedPhi, '\0', kMaxBJets * sizeof(float));
00073 std::memset(ohBJetIPL25Tag, '\0', kMaxBJets * sizeof(float));
00074 std::memset(ohBJetIPL3Tag, '\0', kMaxBJets * sizeof(float));
00075 std::memset(ohBJetIPLooseL25Tag, '\0', kMaxBJets * sizeof(float));
00076 std::memset(ohBJetIPLooseL3Tag, '\0', kMaxBJets * sizeof(float));
00077 std::memset(ohBJetMuL25Tag, '\0', kMaxBJets * sizeof(int));
00078 std::memset(ohBJetMuL3Tag, '\0', kMaxBJets * sizeof(float));
00079 std::memset(ohBJetPerfL25Tag, '\0', kMaxBJets * sizeof(int));
00080 std::memset(ohBJetPerfL3Tag, '\0', kMaxBJets * sizeof(int));
00081 }
00082
00083 HLTBJet::~HLTBJet()
00084 { }
00085
00086 void HLTBJet::setup(const edm::ParameterSet & config, TTree * tree)
00087 {
00088
00089 if (tree) {
00090 tree->Branch("NohBJetL2", & NohBJetL2, "NohBJetL2/I");
00091 tree->Branch("ohBJetL2Energy", ohBJetL2Energy, "ohBJetL2Energy[NohBJetL2]/F");
00092 tree->Branch("ohBJetL2Et", ohBJetL2Et, "ohBJetL2Et[NohBJetL2]/F");
00093 tree->Branch("ohBJetL2Pt", ohBJetL2Pt, "ohBJetL2Pt[NohBJetL2]/F");
00094 tree->Branch("ohBJetL2Eta", ohBJetL2Eta, "ohBJetL2Eta[NohBJetL2]/F");
00095 tree->Branch("ohBJetL2Phi", ohBJetL2Phi, "ohBJetL2Phi[NohBJetL2]/F");
00096
00097 tree->Branch("NohBJetL2Corrected", & NohBJetL2Corrected, "NohBJetL2Corrected/I");
00098 tree->Branch("ohBJetL2CorrectedEnergy", ohBJetL2CorrectedEnergy, "ohBJetL2CorrectedEnergy[NohBJetL2Corrected]/F");
00099 tree->Branch("ohBJetL2CorrectedEt", ohBJetL2CorrectedEt, "ohBJetL2CorrectedEt[NohBJetL2Corrected]/F");
00100 tree->Branch("ohBJetL2CorrectedPt", ohBJetL2CorrectedPt, "ohBJetL2CorrectedPt[NohBJetL2Corrected]/F");
00101 tree->Branch("ohBJetL2CorrectedEta", ohBJetL2CorrectedEta, "ohBJetL2CorrectedEta[NohBJetL2Corrected]/F");
00102 tree->Branch("ohBJetL2CorrectedPhi", ohBJetL2CorrectedPhi, "ohBJetL2CorrectedPhi[NohBJetL2Corrected]/F");
00103
00104 tree->Branch("ohBJetIPL25Tag", ohBJetIPL25Tag, "ohBJetIPL25Tag[NohBJetL2]/F");
00105 tree->Branch("ohBJetIPL3Tag", ohBJetIPL3Tag, "ohBJetIPL3Tag[NohBJetL2]/F");
00106 tree->Branch("ohBJetIPLooseL25Tag", ohBJetIPLooseL25Tag, "ohBJetIPLooseL25Tag[NohBJetL2]/F");
00107 tree->Branch("ohBJetIPLooseL3Tag", ohBJetIPLooseL3Tag, "ohBJetIPLooseL3Tag[NohBJetL2]/F");
00108 tree->Branch("ohBJetMuL25Tag", ohBJetMuL25Tag, "ohBJetMuL25Tag[NohBJetL2]/I");
00109 tree->Branch("ohBJetMuL3Tag", ohBJetMuL3Tag, "ohBJetMuL3Tag[NohBJetL2]/F");
00110 tree->Branch("ohBJetPerfL25Tag", ohBJetPerfL25Tag, "ohBJetPerfL25Tag[NohBJetL2]/I");
00111 tree->Branch("ohBJetPerfL3Tag", ohBJetPerfL3Tag, "ohBJetPerfL3Tag[NohBJetL2]/I");
00112 }
00113 }
00114
00115 void HLTBJet::analyze(
00116 const edm::Handle<edm::View<reco::Jet> > & rawBJets,
00117 const edm::Handle<edm::View<reco::Jet> > & correctedBJets,
00118 const edm::Handle<reco::JetTagCollection> & lifetimeBJetsL25,
00119 const edm::Handle<reco::JetTagCollection> & lifetimeBJetsL3,
00120 const edm::Handle<reco::JetTagCollection> & lifetimeBJetsL25Relaxed,
00121 const edm::Handle<reco::JetTagCollection> & lifetimeBJetsL3Relaxed,
00122 const edm::Handle<reco::JetTagCollection> & softmuonBJetsL25,
00123 const edm::Handle<reco::JetTagCollection> & softmuonBJetsL3,
00124 const edm::Handle<reco::JetTagCollection> & performanceBJetsL25,
00125 const edm::Handle<reco::JetTagCollection> & performanceBJetsL3,
00126 TTree * tree)
00127 {
00128
00129 clear();
00130
00131
00132 if (rawBJets.isValid())
00133 analyseJets(* rawBJets);
00134
00135 if (correctedBJets.isValid())
00136 analyseCorrectedJets(* correctedBJets);
00137
00138 if (rawBJets.isValid() and lifetimeBJetsL25.isValid() and lifetimeBJetsL3.isValid())
00139 analyseLifetime(* rawBJets, * lifetimeBJetsL25, * lifetimeBJetsL3);
00140
00141 if (rawBJets.isValid() and lifetimeBJetsL25Relaxed.isValid() and lifetimeBJetsL3Relaxed.isValid())
00142 analyseLifetimeLoose(* rawBJets, * lifetimeBJetsL25Relaxed, * lifetimeBJetsL3Relaxed);
00143
00144 if (rawBJets.isValid() and softmuonBJetsL25.isValid() and softmuonBJetsL3.isValid())
00145 analyseSoftmuon(* rawBJets, * softmuonBJetsL25, * softmuonBJetsL3);
00146
00147 if (rawBJets.isValid() and performanceBJetsL25.isValid() and performanceBJetsL3.isValid())
00148 analysePerformance(* rawBJets, * performanceBJetsL25, * performanceBJetsL3);
00149 }
00150
00151 void HLTBJet::analyseJets(const edm::View<reco::Jet> & jets)
00152 {
00153 size_t size = std::min(kMaxBJets, jets.size());
00154 NohBJetL2 = size;
00155 for (size_t i = 0; i < size; ++i) {
00156 ohBJetL2Energy[i] = jets[i].energy();
00157 ohBJetL2Et[i] = jets[i].et();
00158 ohBJetL2Pt[i] = jets[i].pt();
00159 ohBJetL2Eta[i] = jets[i].eta();
00160 ohBJetL2Phi[i] = jets[i].phi();
00161 }
00162 }
00163
00164 void HLTBJet::analyseCorrectedJets(const edm::View<reco::Jet> & jets)
00165 {
00166 size_t size = std::min(kMaxBJets, jets.size());
00167 NohBJetL2Corrected = size;
00168 for (size_t i = 0; i < size; ++i) {
00169 ohBJetL2CorrectedEnergy[i] = jets[i].energy();
00170 ohBJetL2CorrectedEt[i] = jets[i].et();
00171 ohBJetL2CorrectedPt[i] = jets[i].pt();
00172 ohBJetL2CorrectedEta[i] = jets[i].eta();
00173 ohBJetL2CorrectedPhi[i] = jets[i].phi();
00174 }
00175 }
00176
00177 void HLTBJet::analyseLifetime(
00178 const edm::View<reco::Jet> & jets,
00179 const reco::JetTagCollection & tagsL25,
00180 const reco::JetTagCollection & tagsL3)
00181 {
00182 if (tagsL25.size() != jets.size()) {
00183 edm::LogWarning("OpenHLT") << kBTagLifetimeBJetsL25 << " collection has " << tagsL25.size() << " elements, but " << jets.size() << " where expected from L2" << std::endl;
00184 return;
00185 }
00186 if (tagsL3.size() != jets.size()) {
00187 edm::LogWarning("OpenHLT") << kBTagLifetimeBJetsL3 << " collection has " << tagsL3.size() << " elements, but " << jets.size() << " where expected from L2" << std::endl;
00188 return;
00189 }
00190 size_t size = std::min(kMaxBJets, jets.size());
00191 for (size_t i = 0; i < size; i++) {
00192 ohBJetIPL25Tag[i] = tagsL25[i].second;
00193 ohBJetIPL3Tag[i] = tagsL3[i].second;
00194 }
00195 }
00196
00197 void HLTBJet::analyseLifetimeLoose(
00198 const edm::View<reco::Jet> & jets,
00199 const reco::JetTagCollection & tagsL25,
00200 const reco::JetTagCollection & tagsL3)
00201 {
00202 if (tagsL25.size() != jets.size()) {
00203 edm::LogWarning("OpenHLT") << kBTagLifetimeBJetsL25Relaxed << " collection has " << tagsL25.size() << " elements, but " << jets.size() << " where expected from L2" << std::endl;
00204 return;
00205 }
00206 if (tagsL3.size() != jets.size()) {
00207 edm::LogWarning("OpenHLT") << kBTagLifetimeBJetsL3Relaxed << " collection has " << tagsL3.size() << " elements, but " << jets.size() << " where expected from L2" << std::endl;
00208 return;
00209 }
00210 size_t size = std::min(kMaxBJets, jets.size());
00211 for (size_t i = 0; i < size; i++) {
00212 ohBJetIPLooseL25Tag[i] = tagsL25[i].second;
00213 ohBJetIPLooseL3Tag[i] = tagsL3[i].second;
00214 }
00215 }
00216
00217 void HLTBJet::analyseSoftmuon(
00218 const edm::View<reco::Jet> & jets,
00219 const reco::JetTagCollection & tagsL25,
00220 const reco::JetTagCollection & tagsL3)
00221 {
00222 if (tagsL25.size() != jets.size()) {
00223 edm::LogWarning("OpenHLT") << kBTagSoftmuonBJetsL25 << " collection has " << tagsL25.size() << " elements, but " << jets.size() << " where expected from L2" << std::endl;
00224 return;
00225 }
00226 if (tagsL3.size() != jets.size()) {
00227 edm::LogWarning("OpenHLT") << kBTagSoftmuonBJetsL3 << " collection has " << tagsL3.size() << " elements, but " << jets.size() << " where expected from L2" << std::endl;
00228 return;
00229 }
00230 size_t size = std::min(kMaxBJets, jets.size());
00231 for (size_t i = 0; i < size; i++) {
00232 ohBJetMuL25Tag[i] = (tagsL25[i].second > 0.) ? 1 : 0;
00233 ohBJetMuL3Tag[i] = tagsL3[i].second;
00234 }
00235 }
00236
00237 void HLTBJet::analysePerformance(
00238 const edm::View<reco::Jet> & jets,
00239 const reco::JetTagCollection & tagsL25,
00240 const reco::JetTagCollection & tagsL3)
00241 {
00242 if (tagsL25.size() != jets.size()) {
00243 edm::LogWarning("OpenHLT") << kBTagPerformanceBJetsL25 << " collection has " << tagsL25.size() << " elements, but " << jets.size() << " where expected from L2" << std::endl;
00244 return;
00245 }
00246 if (tagsL3.size() != jets.size()) {
00247 edm::LogWarning("OpenHLT") << kBTagPerformanceBJetsL3 << " collection has " << tagsL3.size() << " elements, but " << jets.size() << " where expected from L2" << std::endl;
00248 return;
00249 }
00250 size_t size = std::min(kMaxBJets, jets.size());
00251 for (size_t i = 0; i < size; i++) {
00252 ohBJetPerfL25Tag[i] = (tagsL25[i].second > 0.) ? 1 : 0;
00253 ohBJetPerfL3Tag[i] = (tagsL3[i].second > 0.) ? 1 : 0;
00254 }
00255 }