File indexing completed on 2024-04-06 12:32:57
0001 #include "DataFormats/BTauReco/interface/SecondaryVertexTagInfo.h"
0002 #include "TH1D.h"
0003 #include "TH1F.h"
0004 #include "TH2D.h"
0005 #include "TMath.h"
0006 #include <Math/Functions.h>
0007 #include <Math/SMatrix.h>
0008 #include <Math/SVector.h>
0009 #include <cmath>
0010 #include <map>
0011 #include <string>
0012
0013 #include "DQMServices/Core/interface/DQMStore.h"
0014 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0015 #include "FWCore/Framework/interface/Frameworkfwd.h"
0016 #include "FWCore/Framework/interface/MakerMacros.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018 #include "FWCore/ServiceRegistry/interface/Service.h"
0019
0020 #include "CommonTools/Statistics/interface/ChiSquaredProbability.h"
0021 #include "DataFormats/Math/interface/LorentzVector.h"
0022 #include "DataFormats/Math/interface/Vector.h"
0023 #include "Math/VectorUtil.h"
0024 #include "SimTracker/TrackHistory/interface/VertexClassifierByProxy.h"
0025 #include "TROOT.h"
0026 #include <Math/GenVector/PxPyPzE4D.h>
0027 #include <Math/GenVector/PxPyPzM4D.h>
0028 #include <TVector3.h>
0029
0030
0031
0032 using namespace reco;
0033 using namespace std;
0034 using namespace edm;
0035
0036 class recoBSVTagInfoValidationAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0037 public:
0038 typedef dqm::legacy::DQMStore DQMStore;
0039 typedef dqm::legacy::MonitorElement MonitorElement;
0040
0041 explicit recoBSVTagInfoValidationAnalyzer(const edm::ParameterSet &);
0042
0043 private:
0044 void analyze(const edm::Event &, const edm::EventSetup &) override;
0045 void endJob() override;
0046
0047
0048 VertexClassifierByProxy<reco::SecondaryVertexTagInfoCollection> classifier_;
0049
0050 Int_t numberVertexClassifier_;
0051 edm::InputTag trackingTruth_;
0052 edm::InputTag svTagInfoProducer_;
0053
0054 DQMStore *dqmStore_;
0055 std::string dqmLabel;
0056
0057 Int_t n_event;
0058 Int_t rs_total_nall;
0059 Int_t rs_total_nsv;
0060 Int_t rs_total_nbv;
0061 Int_t rs_total_nbsv;
0062 Int_t rs_total_ncv;
0063 Int_t rs_total_nlv;
0064 Int_t total_nfake;
0065
0066 Int_t sr_total_nall;
0067 Int_t sr_total_nsv;
0068 Int_t sr_total_nbv;
0069 Int_t sr_total_nbsv;
0070 Int_t sr_total_ncv;
0071 Int_t sr_total_nlv;
0072 Int_t total_nmiss;
0073
0074
0075 void bookRecoToSim(std::string const &);
0076 void bookSimToReco(std::string const &);
0077
0078
0079 void fillRecoToSim(std::string const &, reco::Vertex const &, TrackingVertexRef const &);
0080 void fillSimToReco(std::string const &, reco::VertexBaseRef const &, TrackingVertexRef const &);
0081
0082
0083 std::map<std::string, MonitorElement *> HistIndex_;
0084
0085
0086 edm::EDGetTokenT<reco::SecondaryVertexTagInfoCollection> svInfoToken;
0087 edm::EDGetTokenT<TrackingVertexCollection> tvToken;
0088 };
0089
0090 recoBSVTagInfoValidationAnalyzer::recoBSVTagInfoValidationAnalyzer(const edm::ParameterSet &config)
0091 : classifier_(config, consumesCollector()) {
0092
0093 n_event = 0;
0094 rs_total_nall = 0;
0095 rs_total_nsv = 0;
0096 rs_total_nbv = 0;
0097 rs_total_nbsv = 0;
0098 rs_total_ncv = 0;
0099 rs_total_nlv = 0;
0100 total_nfake = 0;
0101
0102 sr_total_nall = 0;
0103 sr_total_nsv = 0;
0104 sr_total_nbv = 0;
0105 sr_total_nbsv = 0;
0106 sr_total_ncv = 0;
0107 sr_total_nlv = 0;
0108 total_nmiss = 0;
0109
0110 usesResource("DQMStore");
0111
0112 dqmStore_ = edm::Service<DQMStore>().operator->();
0113 dqmLabel = "SVValidation/";
0114 dqmStore_->setCurrentFolder(dqmLabel);
0115
0116
0117 svInfoToken = consumes<reco::SecondaryVertexTagInfoCollection>(config.getParameter<InputTag>("svTagInfoProducer"));
0118
0119 tvToken = consumes<TrackingVertexCollection>(config.getParameter<InputTag>("trackingTruth"));
0120
0121 numberVertexClassifier_ = VertexCategories::Unknown + 1;
0122
0123
0124 HistIndex_["VertexClassifier"] = dqmStore_->book1D("VertexClassifier",
0125 "Frequency for the different track categories",
0126 numberVertexClassifier_,
0127 -0.5,
0128 numberVertexClassifier_ - 0.5);
0129
0130
0131 HistIndex_["rs_All_MatchQuality"] = dqmStore_->book1D("rs_All_MatchQuality", "Quality of Match", 51, -0.01, 1.01);
0132 HistIndex_["rs_All_FlightDistance2d"] =
0133 dqmStore_->book1D("rs_All_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
0134 HistIndex_["rs_SecondaryVertex_FlightDistance2d"] =
0135 dqmStore_->book1D("rs_SecondaryVertex_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
0136 HistIndex_["rs_BSV_FlightDistance2d"] =
0137 dqmStore_->book1D("rs_BSV_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
0138 HistIndex_["rs_BWeakDecay_FlightDistance2d"] =
0139 dqmStore_->book1D("rs_BWeakDecay_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
0140 HistIndex_["rs_CWeakDecay_FlightDistance2d"] =
0141 dqmStore_->book1D("rs_CWeakDecay_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
0142 HistIndex_["rs_Light_FlightDistance2d"] =
0143 dqmStore_->book1D("rs_Light_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5);
0144
0145 HistIndex_["rs_All_nRecVtx"] = dqmStore_->book1D("rs_All_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0146 HistIndex_["rs_SecondaryVertex_nRecVtx"] =
0147 dqmStore_->book1D("rs_SecondaryVertex_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0148 HistIndex_["rs_BSV_nRecVtx"] = dqmStore_->book1D("rs_BSV_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0149 HistIndex_["rs_BWeakDecay_nRecVtx"] =
0150 dqmStore_->book1D("rs_BWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0151 HistIndex_["rs_CWeakDecay_nRecVtx"] =
0152 dqmStore_->book1D("rs_CWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0153 HistIndex_["rs_Light_nRecVtx"] =
0154 dqmStore_->book1D("rs_Light_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0155
0156
0157 HistIndex_["sr_All_MatchQuality"] = dqmStore_->book1D("sr_All_MatchQuality", "Quality of Match", 51, -0.01, 1.01);
0158 HistIndex_["sr_All_nRecVtx"] = dqmStore_->book1D("sr_All_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0159 HistIndex_["sr_SecondaryVertex_nRecVtx"] =
0160 dqmStore_->book1D("sr_SecondaryVertex_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0161 HistIndex_["sr_BSV_nRecVtx"] = dqmStore_->book1D("sr_BSV_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0162 HistIndex_["sr_BWeakDecay_nRecVtx"] =
0163 dqmStore_->book1D("sr_BWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0164 HistIndex_["sr_CWeakDecay_nRecVtx"] =
0165 dqmStore_->book1D("sr_CWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0166 HistIndex_["sr_Light_nRecVtx"] =
0167 dqmStore_->book1D("sr_Light_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5);
0168
0169
0170 for (Int_t i = 0; i < numberVertexClassifier_; ++i)
0171 HistIndex_["VertexClassifier"]->setBinLabel(i + 1, VertexCategories::Names[i]);
0172
0173
0174 bookRecoToSim("rs_All");
0175 bookRecoToSim("rs_SecondaryVertex");
0176 bookRecoToSim("rs_BSV");
0177 bookRecoToSim("rs_BWeakDecay");
0178 bookRecoToSim("rs_CWeakDecay");
0179 bookRecoToSim("rs_Light");
0180
0181 bookSimToReco("sr_All");
0182 bookSimToReco("sr_SecondaryVertex");
0183 bookSimToReco("sr_BSV");
0184 bookSimToReco("sr_BWeakDecay");
0185 bookSimToReco("sr_CWeakDecay");
0186 bookSimToReco("sr_Light");
0187 }
0188
0189 void recoBSVTagInfoValidationAnalyzer::analyze(const edm::Event &event, const edm::EventSetup &setup) {
0190 ++n_event;
0191
0192 std::cout << "*** Analyzing " << event.id() << " n_event = " << n_event << std::endl << std::endl;
0193
0194
0195 classifier_.newEvent(event, setup);
0196
0197
0198 edm::Handle<reco::SecondaryVertexTagInfoCollection> svTagInfoCollection;
0199 event.getByToken(svInfoToken, svTagInfoCollection);
0200
0201 VertexHistory const &tracer = classifier_.history();
0202
0203 cout << "* Event " << n_event << " ; svTagInfoCollection->size() = " << svTagInfoCollection->size() << endl;
0204
0205 int rs_nall = 0;
0206 int rs_nsv = 0;
0207 int rs_nbv = 0;
0208 int rs_nbsv = 0;
0209 int rs_ncv = 0;
0210 int rs_nlv = 0;
0211 int nfake = 0;
0212
0213 int sr_nall = 0;
0214 int sr_nsv = 0;
0215 int sr_nbv = 0;
0216 int sr_nbsv = 0;
0217 int sr_ncv = 0;
0218 int sr_nlv = 0;
0219 int nmiss = 0;
0220
0221
0222 for (std::size_t index = 0; index < svTagInfoCollection->size(); ++index) {
0223 reco::SecondaryVertexTagInfoRef svTagInfo(svTagInfoCollection, index);
0224
0225
0226 for (std::size_t vindex = 0; vindex < svTagInfo->nVertices(); ++vindex) {
0227
0228 classifier_.evaluate(svTagInfo, vindex);
0229
0230
0231 double rs_quality = tracer.quality();
0232
0233
0234 for (Int_t i = 0; i != numberVertexClassifier_; ++i) {
0235 if (classifier_.is((VertexCategories::Category)i)) {
0236 HistIndex_["VertexClassifier"]->Fill(i);
0237 }
0238 }
0239 if (!classifier_.is(VertexCategories::Fake)) {
0240 HistIndex_["rs_All_MatchQuality"]->Fill(rs_quality);
0241 fillRecoToSim("rs_All", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
0242 HistIndex_["rs_All_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
0243 rs_nall++;
0244
0245 if (classifier_.is(VertexCategories::SecondaryVertex)) {
0246 fillRecoToSim("rs_SecondaryVertex", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
0247 HistIndex_["rs_SecondaryVertex_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
0248 rs_nsv++;
0249 }
0250
0251 if (classifier_.is(VertexCategories::BWeakDecay)) {
0252 fillRecoToSim("rs_BWeakDecay", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
0253 HistIndex_["rs_BWeakDecay_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
0254 rs_nbv++;
0255
0256 if (classifier_.is(VertexCategories::SecondaryVertex)) {
0257 fillRecoToSim("rs_BSV", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
0258 HistIndex_["rs_BSV_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
0259 rs_nbsv++;
0260 }
0261 }
0262
0263 else if (classifier_.is(VertexCategories::CWeakDecay)) {
0264 fillRecoToSim("rs_CWeakDecay", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
0265 HistIndex_["rs_CWeakDecay_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
0266 rs_ncv++;
0267
0268 } else {
0269 fillRecoToSim("rs_Light", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
0270 HistIndex_["rs_Light_FlightDistance2d"]->Fill(svTagInfo->flightDistance(vindex, true).value());
0271 rs_nlv++;
0272 }
0273 }
0274
0275 else {
0276 cout << " VertexCategories::Fake!!" << endl;
0277 nfake++;
0278 }
0279
0280 }
0281
0282 }
0283
0284 HistIndex_["rs_All_nRecVtx"]->Fill(rs_nall);
0285 HistIndex_["rs_SecondaryVertex_nRecVtx"]->Fill(rs_nsv);
0286 HistIndex_["rs_BWeakDecay_nRecVtx"]->Fill(rs_nbv);
0287 HistIndex_["rs_BSV_nRecVtx"]->Fill(rs_nbsv);
0288 HistIndex_["rs_CWeakDecay_nRecVtx"]->Fill(rs_ncv);
0289 HistIndex_["rs_Light_nRecVtx"]->Fill(rs_nlv);
0290 cout << endl;
0291
0292
0293
0294
0295
0296 edm::Handle<TrackingVertexCollection> TVCollection;
0297 event.getByToken(tvToken, TVCollection);
0298
0299 for (std::size_t index = 0; index < TVCollection->size(); ++index) {
0300 TrackingVertexRef trackingVertex(TVCollection, index);
0301
0302 classifier_.evaluate(trackingVertex);
0303
0304 double sr_quality = tracer.quality();
0305
0306 if (classifier_.is(VertexCategories::Reconstructed)) {
0307 HistIndex_["sr_All_MatchQuality"]->Fill(sr_quality);
0308 fillSimToReco("sr_All", tracer.recoVertex(), trackingVertex);
0309 sr_nall++;
0310
0311 if (classifier_.is(VertexCategories::SecondaryVertex)) {
0312 fillSimToReco("sr_SecondaryVertex", tracer.recoVertex(), trackingVertex);
0313 sr_nsv++;
0314 }
0315
0316 if (classifier_.is(VertexCategories::BWeakDecay)) {
0317 fillSimToReco("sr_BWeakDecay", tracer.recoVertex(), trackingVertex);
0318 sr_nbv++;
0319
0320 if (classifier_.is(VertexCategories::SecondaryVertex)) {
0321 fillSimToReco("sr_BSV", tracer.recoVertex(), trackingVertex);
0322 sr_nbsv++;
0323 }
0324
0325 }
0326
0327 else if (classifier_.is(VertexCategories::CWeakDecay)) {
0328 fillSimToReco("sr_CWeakDecay", tracer.recoVertex(), trackingVertex);
0329 sr_ncv++;
0330 }
0331
0332 else {
0333 fillSimToReco("sr_Light", tracer.recoVertex(), trackingVertex);
0334 sr_nlv++;
0335 }
0336
0337 }
0338 else {
0339
0340 nmiss++;
0341 }
0342
0343 }
0344
0345 HistIndex_["sr_All_nRecVtx"]->Fill(sr_nall);
0346 HistIndex_["sr_SecondaryVertex_nRecVtx"]->Fill(sr_nsv);
0347 HistIndex_["sr_BWeakDecay_nRecVtx"]->Fill(sr_nbv);
0348 HistIndex_["sr_BSV_nRecVtx"]->Fill(sr_nbsv);
0349 HistIndex_["sr_CWeakDecay_nRecVtx"]->Fill(sr_ncv);
0350 HistIndex_["sr_Light_nRecVtx"]->Fill(rs_nlv);
0351
0352 rs_total_nall += rs_nall;
0353 rs_total_nsv += rs_nsv;
0354 rs_total_nbv += rs_nbv;
0355 rs_total_nbsv += rs_nbsv;
0356 rs_total_ncv += rs_ncv;
0357 rs_total_nlv += rs_nlv;
0358 total_nfake += nfake;
0359
0360 sr_total_nall += sr_nall;
0361 sr_total_nsv += sr_nsv;
0362 sr_total_nbv += sr_nbv;
0363 sr_total_nbsv += sr_nbsv;
0364 sr_total_ncv += sr_ncv;
0365 sr_total_nlv += sr_nlv;
0366 total_nmiss += nmiss;
0367 }
0368
0369 void recoBSVTagInfoValidationAnalyzer::bookRecoToSim(std::string const &prefix) {
0370
0371
0372 std::string name = prefix + "_Pullx";
0373 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
0374 name = prefix + "_Pully";
0375 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
0376 name = prefix + "_Pullz";
0377 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
0378
0379 name = prefix + "_Resx";
0380 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
0381 name = prefix + "_Resy";
0382 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
0383 name = prefix + "_Resz";
0384 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
0385
0386 name = prefix + "_Chi2Norm";
0387 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, 0, 10.);
0388 name = prefix + "_Chi2Prob";
0389 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, 0., 1.);
0390
0391 name = prefix + "_nRecTrks";
0392 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 501, -0.5, 500.5);
0393
0394 name = prefix + "_AverageTrackWeight";
0395 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.1, 1.1);
0396
0397 name = prefix + "_Mass";
0398 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 65, 0., 6.5);
0399
0400 name = prefix + "_RecPt";
0401 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
0402
0403 name = prefix + "_RecEta";
0404 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
0405
0406 name = prefix + "_RecCharge";
0407 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 21, -0.5, 20.5);
0408
0409 name = prefix + "_RecTrackPt";
0410 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
0411
0412 name = prefix + "_RecTrackEta";
0413 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
0414
0415 name = prefix + "_nSimTrks";
0416 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 501, -0.5, 500.5);
0417
0418 name = prefix + "_SimPt";
0419 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
0420
0421 name = prefix + "_SimEta";
0422 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
0423
0424 name = prefix + "_SimCharge";
0425 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 21, -0.5, 20.5);
0426
0427 name = prefix + "_SimTrackPt";
0428 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 500, 0., 500.);
0429
0430 name = prefix + "_SimTrackEta";
0431 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
0432 }
0433
0434 void recoBSVTagInfoValidationAnalyzer::bookSimToReco(std::string const &prefix) {
0435
0436
0437 std::string name = prefix + "_Pullx";
0438 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
0439 name = prefix + "_Pully";
0440 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
0441 name = prefix + "_Pullz";
0442 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -10., 10.);
0443
0444 name = prefix + "_Resx";
0445 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
0446 name = prefix + "_Resy";
0447 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
0448 name = prefix + "_Resz";
0449 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.05, 0.05);
0450
0451 name = prefix + "_Chi2Norm";
0452 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, 0, 10.);
0453 name = prefix + "_Chi2Prob";
0454 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, 0., 1.);
0455
0456 name = prefix + "_nRecTrks";
0457 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 501, -0.5, 500.5);
0458
0459 name = prefix + "_AverageTrackWeight";
0460 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 100, -0.1, 1.1);
0461
0462 name = prefix + "_Mass";
0463 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 65, 0., 6.5);
0464
0465 name = prefix + "_RecPt";
0466 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
0467
0468 name = prefix + "_RecEta";
0469 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
0470
0471 name = prefix + "_RecCharge";
0472 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 21, -0.5, 20.5);
0473
0474 name = prefix + "_RecTrackPt";
0475 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
0476
0477 name = prefix + "_RecTrackEta";
0478 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
0479
0480 name = prefix + "_nSimTrks";
0481 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 501, -0.5, 500.5);
0482
0483 name = prefix + "_SimPt";
0484 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 2000, 0., 1000.);
0485
0486 name = prefix + "_SimEta";
0487 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
0488
0489 name = prefix + "_SimCharge";
0490 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 21, -0.5, 20.5);
0491
0492 name = prefix + "_SimTrackPt";
0493 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 500, 0., 500.);
0494
0495 name = prefix + "_SimTrackEta";
0496 HistIndex_[name] = dqmStore_->book1D(name.c_str(), name.c_str(), 200, -3., 3.);
0497 }
0498
0499 void recoBSVTagInfoValidationAnalyzer::fillRecoToSim(std::string const &prefix,
0500 reco::Vertex const &vertex,
0501 TrackingVertexRef const &simVertex) {
0502 double pullx = (vertex.x() - simVertex->position().x()) / vertex.xError();
0503 double pully = (vertex.y() - simVertex->position().y()) / vertex.yError();
0504 double pullz = (vertex.z() - simVertex->position().z()) / vertex.zError();
0505
0506 double resx = vertex.x() - simVertex->position().x();
0507 double resy = vertex.y() - simVertex->position().y();
0508 double resz = vertex.z() - simVertex->position().z();
0509
0510 double chi2norm = vertex.normalizedChi2();
0511 double chi2prob = ChiSquaredProbability(vertex.chi2(), vertex.ndof());
0512
0513 double sum_weight = 0.;
0514 double weight = 0.;
0515 double tracksize = vertex.tracksSize();
0516 math::XYZVector momentum;
0517 math::XYZTLorentzVector sum;
0518 int charge = 0;
0519 double thePiMass = 0.13957;
0520 for (reco::Vertex::trackRef_iterator recDaughter = vertex.tracks_begin(); recDaughter != vertex.tracks_end();
0521 ++recDaughter) {
0522 ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<double>> vec;
0523
0524 vec.SetPx((**recDaughter).px());
0525 vec.SetPy((**recDaughter).py());
0526 vec.SetPz((**recDaughter).pz());
0527 vec.SetM(thePiMass);
0528
0529 sum += vec;
0530
0531 weight = vertex.trackWeight(*recDaughter);
0532 sum_weight += weight;
0533
0534 math::XYZVector p;
0535 p = (**recDaughter).momentum();
0536 momentum += p;
0537
0538 charge += (*recDaughter)->charge();
0539
0540 HistIndex_[prefix + "_RecTrackPt"]->Fill((*recDaughter)->pt());
0541 HistIndex_[prefix + "_RecTrackEta"]->Fill((*recDaughter)->eta());
0542 }
0543
0544
0545
0546 double mass = sum.M();
0547 double pt = sqrt(momentum.Perp2());
0548 double eta = momentum.Eta();
0549
0550 math::XYZVector simmomentum;
0551 int simcharge = 0;
0552 for (TrackingVertex::tp_iterator simDaughter = simVertex->daughterTracks_begin();
0553 simDaughter != simVertex->daughterTracks_end();
0554 ++simDaughter) {
0555 math::XYZVector p;
0556 p = (**simDaughter).momentum();
0557 simmomentum += p;
0558
0559 simcharge += (*simDaughter)->charge();
0560
0561 HistIndex_[prefix + "_SimTrackPt"]->Fill((*simDaughter)->pt());
0562 HistIndex_[prefix + "_SimTrackEta"]->Fill((*simDaughter)->eta());
0563 }
0564
0565 double simpt = sqrt(simmomentum.Perp2());
0566 double simeta = simmomentum.Eta();
0567
0568
0569
0570
0571
0572 HistIndex_[prefix + "_nRecTrks"]->Fill(vertex.tracksSize());
0573 HistIndex_[prefix + "_nSimTrks"]->Fill(simVertex->nDaughterTracks());
0574 HistIndex_[prefix + "_Pullx"]->Fill(pullx);
0575 HistIndex_[prefix + "_Pully"]->Fill(pully);
0576 HistIndex_[prefix + "_Pullz"]->Fill(pullz);
0577 HistIndex_[prefix + "_Resx"]->Fill(resx);
0578 HistIndex_[prefix + "_Resy"]->Fill(resy);
0579 HistIndex_[prefix + "_Resz"]->Fill(resz);
0580 HistIndex_[prefix + "_AverageTrackWeight"]->Fill(sum_weight / tracksize);
0581 HistIndex_[prefix + "_Chi2Norm"]->Fill(chi2norm);
0582 HistIndex_[prefix + "_Chi2Prob"]->Fill(chi2prob);
0583 HistIndex_[prefix + "_RecPt"]->Fill(pt);
0584 HistIndex_[prefix + "_RecEta"]->Fill(eta);
0585 HistIndex_[prefix + "_RecCharge"]->Fill(charge);
0586 HistIndex_[prefix + "_Mass"]->Fill(mass);
0587 HistIndex_[prefix + "_SimPt"]->Fill(simpt);
0588 HistIndex_[prefix + "_SimEta"]->Fill(simeta);
0589 HistIndex_[prefix + "_SimCharge"]->Fill(simcharge);
0590 }
0591
0592 void recoBSVTagInfoValidationAnalyzer::fillSimToReco(std::string const &prefix,
0593 reco::VertexBaseRef const &vertex,
0594 TrackingVertexRef const &simVertex) {
0595 double pullx = (vertex->x() - simVertex->position().x()) / vertex->xError();
0596 double pully = (vertex->y() - simVertex->position().y()) / vertex->yError();
0597 double pullz = (vertex->z() - simVertex->position().z()) / vertex->zError();
0598
0599 double resx = vertex->x() - simVertex->position().x();
0600 double resy = vertex->y() - simVertex->position().y();
0601 double resz = vertex->z() - simVertex->position().z();
0602
0603 double chi2norm = vertex->normalizedChi2();
0604 double chi2prob = ChiSquaredProbability(vertex->chi2(), vertex->ndof());
0605
0606 double sum_weight = 0.;
0607 double weight = 0.;
0608 double tracksize = vertex->tracksSize();
0609 math::XYZVector momentum;
0610 math::XYZTLorentzVector sum;
0611 int charge = 0;
0612 double thePiMass = 0.13957;
0613 for (reco::Vertex::trackRef_iterator recDaughter = vertex->tracks_begin(); recDaughter != vertex->tracks_end();
0614 ++recDaughter) {
0615 ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<double>> vec;
0616
0617 vec.SetPx((**recDaughter).px());
0618 vec.SetPy((**recDaughter).py());
0619 vec.SetPz((**recDaughter).pz());
0620 vec.SetM(thePiMass);
0621
0622 sum += vec;
0623
0624 weight = vertex->trackWeight(*recDaughter);
0625 sum_weight += weight;
0626
0627 math::XYZVector p;
0628 p = (**recDaughter).momentum();
0629 momentum += p;
0630
0631 charge += (*recDaughter)->charge();
0632
0633 HistIndex_[prefix + "_RecTrackPt"]->Fill((*recDaughter)->pt());
0634 HistIndex_[prefix + "_RecTrackEta"]->Fill((*recDaughter)->eta());
0635 }
0636
0637
0638
0639 double mass = sum.M();
0640 double pt = sqrt(momentum.Perp2());
0641 double eta = momentum.Eta();
0642
0643 math::XYZVector simmomentum;
0644 int simcharge = 0;
0645 for (TrackingVertex::tp_iterator simDaughter = simVertex->daughterTracks_begin();
0646 simDaughter != simVertex->daughterTracks_end();
0647 ++simDaughter) {
0648 math::XYZVector p;
0649 p = (**simDaughter).momentum();
0650 simmomentum += p;
0651
0652 simcharge += (*simDaughter)->charge();
0653
0654 HistIndex_[prefix + "_SimTrackPt"]->Fill((*simDaughter)->pt());
0655 HistIndex_[prefix + "_SimTrackEta"]->Fill((*simDaughter)->eta());
0656 }
0657
0658 double simpt = sqrt(simmomentum.Perp2());
0659 double simeta = simmomentum.Eta();
0660
0661
0662
0663
0664
0665 HistIndex_[prefix + "_nRecTrks"]->Fill(vertex->tracksSize());
0666 HistIndex_[prefix + "_nSimTrks"]->Fill(simVertex->nDaughterTracks());
0667 HistIndex_[prefix + "_Pullx"]->Fill(pullx);
0668 HistIndex_[prefix + "_Pully"]->Fill(pully);
0669 HistIndex_[prefix + "_Pullz"]->Fill(pullz);
0670 HistIndex_[prefix + "_Resx"]->Fill(resx);
0671 HistIndex_[prefix + "_Resy"]->Fill(resy);
0672 HistIndex_[prefix + "_Resz"]->Fill(resz);
0673 HistIndex_[prefix + "_AverageTrackWeight"]->Fill(sum_weight / tracksize);
0674 HistIndex_[prefix + "_Chi2Norm"]->Fill(chi2norm);
0675 HistIndex_[prefix + "_Chi2Prob"]->Fill(chi2prob);
0676 HistIndex_[prefix + "_RecPt"]->Fill(pt);
0677 HistIndex_[prefix + "_RecEta"]->Fill(eta);
0678 HistIndex_[prefix + "_RecCharge"]->Fill(charge);
0679 HistIndex_[prefix + "_Mass"]->Fill(mass);
0680 HistIndex_[prefix + "_SimPt"]->Fill(simpt);
0681 HistIndex_[prefix + "_SimEta"]->Fill(simeta);
0682 HistIndex_[prefix + "_SimCharge"]->Fill(simcharge);
0683 }
0684
0685 void recoBSVTagInfoValidationAnalyzer::endJob() {
0686 std::cout << std::endl;
0687 std::cout << " ====== Total Number of analyzed events: " << n_event << " ====== " << std::endl;
0688 std::cout << " ====== Total Number of R2S All: " << rs_total_nall << " ====== " << std::endl;
0689 std::cout << " ====== Total Number of R2S SecondaryVertex: " << rs_total_nsv << " ====== " << std::endl;
0690 std::cout << " ====== Total Number of R2S BWeakDecay: " << rs_total_nbv << " ====== " << std::endl;
0691 std::cout << " ====== Total Number of R2S BWeakDecay::SecondaryVertex: " << rs_total_nbsv << " ====== " << std::endl;
0692 std::cout << " ====== Total Number of R2S CWeakDecay: " << rs_total_ncv << " ====== " << std::endl;
0693 std::cout << " ====== Total Number of R2S Light: " << rs_total_nlv << " ====== " << std::endl;
0694 std::cout << std::endl;
0695 std::cout << " ====== Total Number of S2R All: " << sr_total_nall << " ====== " << std::endl;
0696 std::cout << " ====== Total Number of S2R SecondaryVertex: " << sr_total_nsv << " ====== " << std::endl;
0697 std::cout << " ====== Total Number of S2R BWeakDecay: " << sr_total_nbv << " ====== " << std::endl;
0698 std::cout << " ====== Total Number of S2R BWeakDecay::SecondaryVertex: " << sr_total_nbsv << " ====== " << std::endl;
0699 std::cout << " ====== Total Number of S2R CWeakDecay: " << sr_total_ncv << " ====== " << std::endl;
0700 std::cout << " ====== Total Number of S2R Light: " << sr_total_nlv << " ====== " << std::endl;
0701 std::cout << std::endl;
0702 std::cout << " ====== Total Number of Fake Vertices: " << total_nfake << " ====== " << std::endl;
0703 }
0704
0705 DEFINE_FWK_MODULE(recoBSVTagInfoValidationAnalyzer);