Back to home page

Project CMSSW displayed by LXR
 
 

    

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 // class decleration
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   // Member data
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   // Bookeeping of all the histograms per category
0075   void bookRecoToSim(std::string const &);
0076   void bookSimToReco(std::string const &);
0077 
0078   // Fill all histogram per category
0079   void fillRecoToSim(std::string const &, reco::Vertex const &, TrackingVertexRef const &);
0080   void fillSimToReco(std::string const &, reco::VertexBaseRef const &, TrackingVertexRef const &);
0081 
0082   // Histogram handlers
0083   std::map<std::string, MonitorElement *> HistIndex_;
0084 
0085   // consumes
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   // Initialize counters
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   //  get the store
0112   dqmStore_ = edm::Service<DQMStore>().operator->();
0113   dqmLabel = "SVValidation/";
0114   dqmStore_->setCurrentFolder(dqmLabel);
0115 
0116   // Get the track collection
0117   svInfoToken = consumes<reco::SecondaryVertexTagInfoCollection>(config.getParameter<InputTag>("svTagInfoProducer"));
0118   // Name of the traking pariticle collection
0119   tvToken = consumes<TrackingVertexCollection>(config.getParameter<InputTag>("trackingTruth"));
0120   // Number of track categories
0121   numberVertexClassifier_ = VertexCategories::Unknown + 1;
0122 
0123   // Define histogram for counting categories
0124   HistIndex_["VertexClassifier"] = dqmStore_->book1D("VertexClassifier",
0125                                                      "Frequency for the different track categories",
0126                                                      numberVertexClassifier_,
0127                                                      -0.5,
0128                                                      numberVertexClassifier_ - 0.5);
0129 
0130   //--- RecoToSim
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   //--- SimToReco
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   // Set the proper categories names
0170   for (Int_t i = 0; i < numberVertexClassifier_; ++i)
0171     HistIndex_["VertexClassifier"]->setBinLabel(i + 1, VertexCategories::Names[i]);
0172 
0173   // book histograms
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   // Set the classifier for a new event
0195   classifier_.newEvent(event, setup);
0196 
0197   // Vertex collection
0198   edm::Handle<reco::SecondaryVertexTagInfoCollection> svTagInfoCollection;
0199   event.getByToken(svInfoToken, svTagInfoCollection);
0200   // Get a constant reference to the track history associated to the classifier
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   // Loop over the svTagInfo collection.
0222   for (std::size_t index = 0; index < svTagInfoCollection->size(); ++index) {
0223     reco::SecondaryVertexTagInfoRef svTagInfo(svTagInfoCollection, index);
0224 
0225     // Loop over the vertexes in svTagInfo
0226     for (std::size_t vindex = 0; vindex < svTagInfo->nVertices(); ++vindex) {
0227       // Classify the vertices
0228       classifier_.evaluate(svTagInfo, vindex);
0229 
0230       // quality of the match
0231       double rs_quality = tracer.quality();
0232 
0233       // Fill the histogram with the categories
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         }  // BWeakDecay
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       }  // end if classifier
0274 
0275       else {
0276         cout << "    VertexCategories::Fake!!" << endl;
0277         nfake++;
0278       }
0279 
0280     }  // end loop over vertices in svTagInfo
0281 
0282   }  // loop over svTagInfo
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   // SIM TO RECO!
0294 
0295   // Vertex collection
0296   edm::Handle<TrackingVertexCollection> TVCollection;
0297   event.getByToken(tvToken, TVCollection);
0298   // Loop over the TV collection.
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       }  // BWeakDecay
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     }  // Reconstructed
0338     else {
0339       // cout << "##### Not reconstructed!" << endl;
0340       nmiss++;
0341     }
0342 
0343   }  // TVCollection.size()
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   // Book pull histograms
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   // Book pull histograms
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   }  // end loop to recDaughters
0543   // cout << "                   average sum of weights = " <<
0544   // sum_weight/tracksize << endl;
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   // cout << "[fillRecoToSim]  vertex.tracksSize() = " << vertex.tracksSize() <<
0569   // " ; simVertex->nDaughterTracks() = " << simVertex->nDaughterTracks() <<
0570   // endl;
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   // cout << "                   average sum of weights = " <<
0637   // sum_weight/tracksize << endl;
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   // cout << "[fillSimToReco]  vertex->tracksSize() = " << vertex->tracksSize()
0662   // << " ; simVertex->nDaughterTracks() = " << simVertex->nDaughterTracks() <<
0663   // endl;
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);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!