Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:43:56

0001 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0002 #include "DQM/TrackingMonitorSource/interface/TrackToTrackComparisonHists.h"
0003 
0004 #include "DQMServices/Core/interface/DQMStore.h"
0005 
0006 #include "DataFormats/TrackReco/interface/Track.h"
0007 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0008 #include "DataFormats/VertexReco/interface/Vertex.h"
0009 #include "DataFormats/Luminosity/interface/LumiSummary.h"
0010 
0011 #include "DataFormats/Math/interface/deltaR.h"
0012 #include "CommonTools/TriggerUtils/interface/GenericTriggerEventFlag.h"
0013 
0014 //
0015 // constructors and destructor
0016 //
0017 TrackToTrackComparisonHists::TrackToTrackComparisonHists(const edm::ParameterSet& iConfig)
0018     : monitoredTrackInputTag_(iConfig.getParameter<edm::InputTag>("monitoredTrack")),
0019       referenceTrackInputTag_(iConfig.getParameter<edm::InputTag>("referenceTrack")),
0020       topDirName_(iConfig.getParameter<std::string>("topDirName")),
0021       dRmin_(iConfig.getParameter<double>("dRmin")),
0022       pTCutForPlateau_(iConfig.getParameter<double>("pTCutForPlateau")),
0023       dxyCutForPlateau_(iConfig.getParameter<double>("dxyCutForPlateau")),
0024       dzWRTPvCut_(iConfig.getParameter<double>("dzWRTPvCut")),
0025       requireValidHLTPaths_(iConfig.getParameter<bool>("requireValidHLTPaths")),
0026       genTriggerEventFlag_(new GenericTriggerEventFlag(
0027           iConfig.getParameter<edm::ParameterSet>("genericTriggerEventPSet"), consumesCollector(), *this))
0028 
0029 {
0030   initialize_parameter(iConfig);
0031 
0032   //now do what ever initialization is needed
0033   monitoredTrackToken_ = consumes<reco::TrackCollection>(monitoredTrackInputTag_);
0034   referenceTrackToken_ = consumes<reco::TrackCollection>(referenceTrackInputTag_);
0035   monitoredBSToken_ = consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("monitoredBeamSpot"));
0036   referenceBSToken_ = consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("referenceBeamSpot"));
0037   monitoredPVToken_ = consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("monitoredPrimaryVertices"));
0038   referencePVToken_ = consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("referencePrimaryVertices"));
0039   lumiScalersToken_ = consumes<LumiScalersCollection>(iConfig.getParameter<edm::InputTag>("scalers"));
0040   onlineMetaDataDigisToken_ =
0041       consumes<OnlineLuminosityRecord>(iConfig.getParameter<edm::InputTag>("onlineMetaDataDigis"));
0042 
0043   referenceTracksMEs_.label = referenceTrackInputTag_.label();
0044   matchedReferenceTracksMEs_.label = referenceTrackInputTag_.label() + "_matched";
0045 
0046   monitoredTracksMEs_.label = monitoredTrackInputTag_.label();
0047   unMatchedMonitoredTracksMEs_.label = monitoredTrackInputTag_.label() + "_unMatched";
0048 
0049   matchTracksMEs_.label = "matches";
0050 }
0051 
0052 TrackToTrackComparisonHists::~TrackToTrackComparisonHists() {
0053   if (genTriggerEventFlag_)
0054     genTriggerEventFlag_.reset();
0055 }
0056 
0057 void TrackToTrackComparisonHists::beginJob(const edm::EventSetup& iSetup) {}
0058 
0059 void TrackToTrackComparisonHists::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0060   LogDebug("TrackToTrackComparisonHists")
0061       << " requireValidHLTPaths_ " << requireValidHLTPaths_ << " hltPathsAreValid_  " << hltPathsAreValid_ << "\n";
0062   // if valid HLT paths are required,
0063   // analyze event only if paths are valid
0064   if (requireValidHLTPaths_ and (not hltPathsAreValid_)) {
0065     return;
0066   }
0067 
0068   LogDebug("TrackToTrackComparisonHists") << " genTriggerEventFlag_->on() " << genTriggerEventFlag_->on()
0069                                           << "  accept:  " << genTriggerEventFlag_->accept(iEvent, iSetup) << "\n";
0070   // Filter out events if Trigger Filtering is requested
0071   if (genTriggerEventFlag_->on() && !genTriggerEventFlag_->accept(iEvent, iSetup)) {
0072     return;
0073   }
0074 
0075   //
0076   //  Get Lumi/LS Info
0077   //
0078 
0079   unsigned int ls = iEvent.id().luminosityBlock();
0080 
0081   double onlinelumi = -1.f;
0082   double PU = -1.f;
0083 
0084   auto const lumiScalersHandle = iEvent.getHandle(lumiScalersToken_);
0085   auto const onlineMetaDataDigisHandle = iEvent.getHandle(onlineMetaDataDigisToken_);
0086 
0087   if (onlineMetaDataDigisHandle.isValid()) {
0088     onlinelumi = onlineMetaDataDigisHandle->instLumi();
0089     PU = onlineMetaDataDigisHandle->avgPileUp();
0090   } else if (lumiScalersHandle.isValid() and not lumiScalersHandle->empty()) {
0091     edm::LogError("TrackToTrackComparisonHists") << "onlineMetaDataDigisHandle not found, trying SCAL";
0092     auto const scalit = lumiScalersHandle->begin();
0093     onlinelumi = scalit->instantLumi();
0094     PU = scalit->pileup();
0095   } else {
0096     edm::LogError("TrackToTrackComparisonHists") << "lumiScalersHandle not found or empty, skipping event";
0097     return;
0098   }
0099 
0100   //
0101   //  Get Reference Track Info
0102   //
0103   edm::Handle<reco::TrackCollection> referenceTracksHandle;
0104   iEvent.getByToken(referenceTrackToken_, referenceTracksHandle);
0105   if (!referenceTracksHandle.isValid()) {
0106     edm::LogError("TrackToTrackComparisonHists") << "referenceTracksHandle not found, skipping event";
0107     return;
0108   }
0109   reco::TrackCollection referenceTracks = *referenceTracksHandle;
0110 
0111   edm::Handle<reco::BeamSpot> referenceBSHandle;
0112   iEvent.getByToken(referenceBSToken_, referenceBSHandle);
0113   if (!referenceBSHandle.isValid()) {
0114     edm::LogError("TrackToTrackComparisonHists") << "referenceBSHandle not found, skipping event";
0115     return;
0116   }
0117   reco::BeamSpot referenceBS = *referenceBSHandle;
0118 
0119   edm::Handle<reco::VertexCollection> referencePVHandle;
0120   iEvent.getByToken(referencePVToken_, referencePVHandle);
0121   if (!referencePVHandle.isValid()) {
0122     edm::LogError("TrackToTrackComparisonHists") << "referencePVHandle not found, skipping event";
0123     return;
0124   }
0125   if (referencePVHandle->empty()) {
0126     edm::LogInfo("TrackToTrackComparisonHists") << "referencePVHandle->size is 0 ";
0127     return;
0128   }
0129   reco::Vertex referencePV = referencePVHandle->at(0);
0130 
0131   //
0132   //  Get Monitored Track Info
0133   //
0134   edm::Handle<reco::TrackCollection> monitoredTracksHandle;
0135   iEvent.getByToken(monitoredTrackToken_, monitoredTracksHandle);
0136   if (!monitoredTracksHandle.isValid()) {
0137     edm::LogError("TrackToTrackComparisonHists") << "monitoredTracksHandle not found, skipping event";
0138     return;
0139   }
0140   reco::TrackCollection monitoredTracks = *monitoredTracksHandle;
0141 
0142   edm::Handle<reco::BeamSpot> monitoredBSHandle;
0143   iEvent.getByToken(monitoredBSToken_, monitoredBSHandle);
0144   if (!monitoredTracksHandle.isValid()) {
0145     edm::LogError("TrackToTrackComparisonHists") << "monitoredBSHandle not found, skipping event";
0146     return;
0147   }
0148   reco::BeamSpot monitoredBS = *monitoredBSHandle;
0149 
0150   edm::Handle<reco::VertexCollection> monitoredPVHandle;
0151   iEvent.getByToken(monitoredPVToken_, monitoredPVHandle);
0152   if (!monitoredPVHandle.isValid()) {
0153     edm::LogError("TrackToTrackComparisonHists") << "monitoredPVHandle not found, skipping event";
0154     return;
0155   }
0156   if (monitoredPVHandle->empty()) {
0157     edm::LogInfo("TrackToTrackComparisonHists") << "monitoredPVHandle->size is 0 ";
0158     return;
0159   }
0160   reco::Vertex monitoredPV = monitoredPVHandle->at(0);
0161 
0162   edm::LogInfo("TrackToTrackComparisonHists")
0163       << "analyzing " << monitoredTrackInputTag_.process() << ":" << monitoredTrackInputTag_.label() << ":"
0164       << monitoredTrackInputTag_.instance() << " w.r.t. " << referenceTrackInputTag_.process() << ":"
0165       << referenceTrackInputTag_.label() << ":" << referenceTrackInputTag_.instance() << " \n";
0166 
0167   //
0168   // Build the dR maps
0169   //
0170   idx2idxByDoubleColl monitored2referenceColl;
0171   fillMap(monitoredTracks, referenceTracks, monitored2referenceColl, dRmin_);
0172 
0173   idx2idxByDoubleColl reference2monitoredColl;
0174   fillMap(referenceTracks, monitoredTracks, reference2monitoredColl, dRmin_);
0175 
0176   unsigned int nReferenceTracks(0);           // Counts the number of refernce tracks
0177   unsigned int nMatchedReferenceTracks(0);    // Counts the number of matched refernce tracks
0178   unsigned int nMonitoredTracks(0);           // Counts the number of monitored tracks
0179   unsigned int nUnmatchedMonitoredTracks(0);  // Counts the number of unmatched monitored tracks
0180 
0181   //
0182   // loop over reference tracks
0183   //
0184   LogDebug("TrackToTrackComparisonHists") << "\n# of tracks (reference): " << referenceTracks.size() << "\n";
0185   for (idx2idxByDoubleColl::const_iterator pItr = reference2monitoredColl.begin(), eItr = reference2monitoredColl.end();
0186        pItr != eItr;
0187        ++pItr) {
0188     nReferenceTracks++;
0189     int trackIdx = pItr->first;
0190     reco::Track track = referenceTracks.at(trackIdx);
0191 
0192     float dzWRTpv = track.dz(referencePV.position());
0193     if (fabs(dzWRTpv) > dzWRTPvCut_)
0194       continue;
0195 
0196     fill_generic_tracks_histos(*&referenceTracksMEs_, &track, &referenceBS, &referencePV, ls, onlinelumi, PU);
0197 
0198     std::map<double, int> trackDRmap = pItr->second;
0199     if (trackDRmap.empty()) {
0200       (matchedReferenceTracksMEs_.h_dRmin)->Fill(-1.);
0201       (matchedReferenceTracksMEs_.h_dRmin_l)->Fill(-1.);
0202       continue;
0203     }
0204 
0205     double dRmin = trackDRmap.begin()->first;
0206     (referenceTracksMEs_.h_dRmin)->Fill(dRmin);
0207     (referenceTracksMEs_.h_dRmin_l)->Fill(dRmin);
0208 
0209     bool matched = false;
0210     if (dRmin < dRmin_)
0211       matched = true;
0212 
0213     if (matched) {
0214       nMatchedReferenceTracks++;
0215       fill_generic_tracks_histos(*&matchedReferenceTracksMEs_, &track, &referenceBS, &referencePV, ls, onlinelumi, PU);
0216       (matchedReferenceTracksMEs_.h_dRmin)->Fill(dRmin);
0217       (matchedReferenceTracksMEs_.h_dRmin_l)->Fill(dRmin);
0218 
0219       int matchedTrackIndex = trackDRmap[dRmin];
0220       reco::Track matchedTrack = monitoredTracks.at(matchedTrackIndex);
0221       fill_matching_tracks_histos(*&matchTracksMEs_, &track, &matchedTrack, &referenceBS, &referencePV);
0222     }
0223 
0224   }  // Over reference tracks
0225 
0226   //
0227   // loop over monitoed tracks
0228   //
0229   LogDebug("TrackToTrackComparisonHists") << "\n# of tracks (monitored): " << monitoredTracks.size() << "\n";
0230   for (idx2idxByDoubleColl::const_iterator pItr = monitored2referenceColl.begin(), eItr = monitored2referenceColl.end();
0231        pItr != eItr;
0232        ++pItr) {
0233     nMonitoredTracks++;
0234     int trackIdx = pItr->first;
0235     reco::Track track = monitoredTracks.at(trackIdx);
0236 
0237     float dzWRTpv = track.dz(monitoredPV.position());
0238     if (fabs(dzWRTpv) > dzWRTPvCut_)
0239       continue;
0240 
0241     fill_generic_tracks_histos(*&monitoredTracksMEs_, &track, &monitoredBS, &monitoredPV, ls, onlinelumi, PU);
0242 
0243     std::map<double, int> trackDRmap = pItr->second;
0244     if (trackDRmap.empty()) {
0245       (unMatchedMonitoredTracksMEs_.h_dRmin)->Fill(-1.);
0246       (unMatchedMonitoredTracksMEs_.h_dRmin_l)->Fill(-1.);
0247       continue;
0248     }
0249 
0250     double dRmin = trackDRmap.begin()->first;
0251     (monitoredTracksMEs_.h_dRmin)->Fill(dRmin);
0252     (monitoredTracksMEs_.h_dRmin_l)->Fill(dRmin);
0253 
0254     bool matched = false;
0255     if (dRmin < dRmin_)
0256       matched = true;
0257 
0258     if (!matched) {
0259       nUnmatchedMonitoredTracks++;
0260       fill_generic_tracks_histos(
0261           *&unMatchedMonitoredTracksMEs_, &track, &monitoredBS, &monitoredPV, ls, onlinelumi, PU);
0262       (unMatchedMonitoredTracksMEs_.h_dRmin)->Fill(dRmin);
0263       (unMatchedMonitoredTracksMEs_.h_dRmin_l)->Fill(dRmin);
0264     }
0265 
0266   }  // over monitoed tracks
0267 
0268   edm::LogInfo("TrackToTrackComparisonHists")
0269       << "Total reference tracks: " << nReferenceTracks << "\n"
0270       << "Total matched reference tracks: " << nMatchedReferenceTracks << "\n"
0271       << "Total monitored tracks: " << nMonitoredTracks << "\n"
0272       << "Total unMatched monitored tracks: " << nUnmatchedMonitoredTracks << "\n";
0273 }
0274 
0275 void TrackToTrackComparisonHists::bookHistograms(DQMStore::IBooker& ibooker,
0276                                                  edm::Run const& iRun,
0277                                                  edm::EventSetup const& iSetup) {
0278   if (genTriggerEventFlag_ && genTriggerEventFlag_->on())
0279     genTriggerEventFlag_->initRun(iRun, iSetup);
0280 
0281   // check if every HLT path specified has a valid match in the HLT Menu
0282   hltPathsAreValid_ =
0283       (genTriggerEventFlag_ && genTriggerEventFlag_->on() && genTriggerEventFlag_->allHLTPathsAreValid());
0284 
0285   // if valid HLT paths are required,
0286   // create DQM outputs only if all paths are valid
0287   if (requireValidHLTPaths_ and (not hltPathsAreValid_)) {
0288     return;
0289   }
0290 
0291   std::string dir = topDirName_;
0292 
0293   bookHistos(ibooker, referenceTracksMEs_, "ref", dir);
0294   bookHistos(ibooker, matchedReferenceTracksMEs_, "ref_matched", dir);
0295 
0296   bookHistos(ibooker, monitoredTracksMEs_, "mon", dir);
0297   bookHistos(ibooker, unMatchedMonitoredTracksMEs_, "mon_unMatched", dir);
0298 
0299   book_matching_tracks_histos(ibooker, matchTracksMEs_, "matches", dir);
0300 }
0301 
0302 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0303 void TrackToTrackComparisonHists::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0304   edm::ParameterSetDescription desc;
0305 
0306   desc.add<bool>("requireValidHLTPaths", true);
0307 
0308   desc.add<edm::InputTag>("monitoredTrack", edm::InputTag("hltMergedTracks"));
0309   desc.add<edm::InputTag>("monitoredBeamSpot", edm::InputTag("hltOnlineBeamSpot"));
0310   desc.add<edm::InputTag>("monitoredPrimaryVertices", edm::InputTag("hltVerticesPFSelector"));
0311 
0312   desc.add<edm::InputTag>("referenceTrack", edm::InputTag("generalTracks"));
0313   desc.add<edm::InputTag>("referenceBeamSpot", edm::InputTag("offlineBeamSpot"));
0314   desc.add<edm::InputTag>("referencePrimaryVertices", edm::InputTag("offlinePrimaryVertices"));
0315   desc.add<edm::InputTag>("scalers", edm::InputTag("scalersRawToDigi"));
0316   desc.add<edm::InputTag>("onlineMetaDataDigis", edm::InputTag("onlineMetaDataDigis"));
0317 
0318   desc.add<std::string>("topDirName", "HLT/Tracking/ValidationWRTOffline");
0319   desc.add<double>("dRmin", 0.002);
0320 
0321   desc.add<double>("pTCutForPlateau", 0.9);
0322   desc.add<double>("dxyCutForPlateau", 2.5);
0323   desc.add<double>("dzWRTPvCut", 1e6);
0324 
0325   edm::ParameterSetDescription genericTriggerEventPSet;
0326   GenericTriggerEventFlag::fillPSetDescription(genericTriggerEventPSet);
0327   desc.add<edm::ParameterSetDescription>("genericTriggerEventPSet", genericTriggerEventPSet);
0328 
0329   edm::ParameterSetDescription histoPSet;
0330   fillHistoPSetDescription(histoPSet);
0331   desc.add<edm::ParameterSetDescription>("histoPSet", histoPSet);
0332 
0333   descriptions.add("trackToTrackComparisonHists", desc);
0334 }
0335 
0336 void TrackToTrackComparisonHists::fillMap(reco::TrackCollection tracks1,
0337                                           reco::TrackCollection tracks2,
0338                                           idx2idxByDoubleColl& map,
0339                                           float dRMin) {
0340   //
0341   // loop on tracks1
0342   //
0343   int i = 0;
0344   for (const auto& track1 : tracks1) {
0345     std::map<double, int> tmp;
0346     int j = 0;
0347     float smallest_dR = 1e9;
0348     int smallest_dR_j = -1;
0349 
0350     //
0351     // loop on tracks2
0352     //
0353     for (const auto& track2 : tracks2) {
0354       double dR = reco::deltaR(track1.eta(), track1.phi(), track2.eta(), track2.phi());
0355 
0356       if (dR < smallest_dR) {
0357         smallest_dR = dR;
0358         smallest_dR_j = j;
0359       }
0360 
0361       if (dR < dRMin) {
0362         tmp[dR] = j;
0363       }
0364 
0365       j++;
0366     }
0367 
0368     //
0369     // If there are no tracks that pass the dR store the smallest (for debugging/validating matching)
0370     //
0371     if (tmp.empty())
0372       tmp[smallest_dR] = smallest_dR_j;
0373 
0374     map.push_back(std::make_pair(i, tmp));
0375     i++;
0376   }
0377 }
0378 
0379 void TrackToTrackComparisonHists::bookHistos(DQMStore::IBooker& ibooker,
0380                                              generalME& mes,
0381                                              TString label,
0382                                              std::string& dir) {
0383   book_generic_tracks_histos(ibooker, mes, label, dir);
0384 }
0385 
0386 void TrackToTrackComparisonHists::book_generic_tracks_histos(DQMStore::IBooker& ibooker,
0387                                                              generalME& mes,
0388                                                              TString label,
0389                                                              std::string& dir) {
0390   ibooker.cd();
0391   ibooker.setCurrentFolder(dir);
0392   (mes.h_pt) = ibooker.book1D(label + "_pt", "track p_{T}", Pt_nbin, Pt_rangeMin, Pt_rangeMax);
0393   (mes.h_eta) = ibooker.book1D(label + "_eta", "track pseudorapidity", Eta_nbin, Eta_rangeMin, Eta_rangeMax);
0394   (mes.h_phi) = ibooker.book1D(label + "_phi", "track #phi", Phi_nbin, Phi_rangeMin, Phi_rangeMax);
0395   (mes.h_dxy) =
0396       ibooker.book1D(label + "_dxy", "track transverse dca to beam spot", Dxy_nbin, Dxy_rangeMin, Dxy_rangeMax);
0397   (mes.h_dz) = ibooker.book1D(label + "_dz", "track longitudinal dca to beam spot", Dz_nbin, Dz_rangeMin, Dz_rangeMax);
0398   (mes.h_dxyWRTpv) = ibooker.book1D(
0399       label + "_dxyWRTpv", "track transverse dca to primary vertex", Dxy_nbin, Dxy_rangeMin, Dxy_rangeMax);
0400   (mes.h_dzWRTpv) = ibooker.book1D(
0401       label + "_dzWRTpv", "track longitudinal dca to primary vertex", Dz_nbin, 0.1 * Dz_rangeMin, 0.1 * Dz_rangeMax);
0402   (mes.h_charge) = ibooker.book1D(label + "_charge", "track charge", 5, -2, 2);
0403   (mes.h_hits) = ibooker.book1D(label + "_hits", "track number of hits", 35, -0.5, 34.5);
0404   (mes.h_dRmin) = ibooker.book1D(label + "_dRmin", "track min dR", 100, 0., 0.01);
0405   (mes.h_dRmin_l) = ibooker.book1D(label + "_dRmin_l", "track min dR", 100, 0., 0.4);
0406 
0407   (mes.h_pt_vs_eta) = ibooker.book2D(label + "_ptVSeta",
0408                                      "track p_{T} vs #eta",
0409                                      Eta_nbin,
0410                                      Eta_rangeMin,
0411                                      Eta_rangeMax,
0412                                      Pt_nbin,
0413                                      Pt_rangeMin,
0414                                      Pt_rangeMax);
0415 
0416   // counts of tracks vs lumi
0417   // for this moment, xmin,xmax and binning are hardcoded, maybe in future in a config file!
0418   // have to add (declare) this in the .h file as well
0419   (mes.h_onlinelumi) = ibooker.book1D(label + "_onlinelumi",
0420                                       "number of tracks vs onlinelumi",
0421                                       onlinelumi_nbin,
0422                                       onlinelumi_rangeMin,
0423                                       onlinelumi_rangeMax);
0424   (mes.h_ls) = ibooker.book1D(label + "_ls", "number of tracks vs ls", ls_nbin, ls_rangeMin, ls_rangeMax);
0425   (mes.h_PU) = ibooker.book1D(label + "_PU", "number of tracks vs PU", PU_nbin, PU_rangeMin, PU_rangeMax);
0426 }
0427 
0428 void TrackToTrackComparisonHists::book_matching_tracks_histos(DQMStore::IBooker& ibooker,
0429                                                               matchingME& mes,
0430                                                               TString label,
0431                                                               std::string& dir) {
0432   ibooker.cd();
0433   ibooker.setCurrentFolder(dir);
0434 
0435   (mes.h_hits_vs_hits) = ibooker.book2D(
0436       label + "_hits_vs_hits", "monitored track # hits vs reference track # hits", 35, -0.5, 34.5, 35, -0.5, 34.5);
0437   (mes.h_pt_vs_pt) = ibooker.book2D(label + "_pt_vs_pt",
0438                                     "monitored track p_{T} vs reference track p_{T}",
0439                                     Pt_nbin,
0440                                     Pt_rangeMin,
0441                                     Pt_rangeMax,
0442                                     Pt_nbin,
0443                                     Pt_rangeMin,
0444                                     Pt_rangeMax);
0445   (mes.h_eta_vs_eta) = ibooker.book2D(label + "_eta_vs_eta",
0446                                       "monitored track #eta vs reference track #eta",
0447                                       Eta_nbin,
0448                                       Eta_rangeMin,
0449                                       Eta_rangeMax,
0450                                       Eta_nbin,
0451                                       Eta_rangeMin,
0452                                       Eta_rangeMax);
0453   (mes.h_phi_vs_phi) = ibooker.book2D(label + "_phi_vs_phi",
0454                                       "monitored track #phi vs reference track #phi",
0455                                       Phi_nbin,
0456                                       Phi_rangeMin,
0457                                       Phi_rangeMax,
0458                                       Phi_nbin,
0459                                       Phi_rangeMin,
0460                                       Phi_rangeMax);
0461 
0462   (mes.h_dPt) = ibooker.book1D(label + "_dPt", "#Delta track #P_T", ptRes_nbin, ptRes_rangeMin, ptRes_rangeMax);
0463   (mes.h_dEta) = ibooker.book1D(label + "_dEta", "#Delta track #eta", etaRes_nbin, etaRes_rangeMin, etaRes_rangeMax);
0464   (mes.h_dPhi) = ibooker.book1D(label + "_dPhi", "#Delta track #phi", phiRes_nbin, phiRes_rangeMin, phiRes_rangeMax);
0465   (mes.h_dDxy) = ibooker.book1D(
0466       label + "_dDxy", "#Delta track transverse dca to beam spot", dxyRes_nbin, dxyRes_rangeMin, dxyRes_rangeMax);
0467   (mes.h_dDz) = ibooker.book1D(
0468       label + "_dDz", "#Delta track longitudinal dca to beam spot", dzRes_nbin, dzRes_rangeMin, dzRes_rangeMax);
0469   (mes.h_dDxyWRTpv) = ibooker.book1D(label + "_dDxyWRTpv",
0470                                      "#Delta track transverse dca to primary vertex ",
0471                                      dxyRes_nbin,
0472                                      dxyRes_rangeMin,
0473                                      dxyRes_rangeMax);
0474   (mes.h_dDzWRTpv) = ibooker.book1D(label + "_dDzWRTpv",
0475                                     "#Delta track longitudinal dca to primary vertex",
0476                                     dzRes_nbin,
0477                                     dzRes_rangeMin,
0478                                     dzRes_rangeMax);
0479   (mes.h_dCharge) = ibooker.book1D(label + "_dCharge", "#Delta track charge", 5, -2.5, 2.5);
0480   (mes.h_dHits) = ibooker.book1D(label + "_dHits", "#Delta track number of hits", 39, -19.5, 19.5);
0481 }
0482 
0483 void TrackToTrackComparisonHists::fill_generic_tracks_histos(generalME& mes,
0484                                                              reco::Track* trk,
0485                                                              reco::BeamSpot* bs,
0486                                                              reco::Vertex* pv,
0487                                                              unsigned int ls,
0488                                                              double onlinelumi,
0489                                                              double PU,
0490                                                              bool requirePlateau) {
0491   float pt = trk->pt();
0492   float eta = trk->eta();
0493   float phi = trk->phi();
0494   float dxy = trk->dxy(bs->position());
0495   float dz = trk->dz(bs->position());
0496   float dxyWRTpv = trk->dxy(pv->position());
0497   float dzWRTpv = trk->dz(pv->position());
0498   float charge = trk->charge();
0499   float nhits = trk->hitPattern().numberOfValidHits();
0500 
0501   bool dxyOnPlateau = (fabs(dxyWRTpv) < dxyCutForPlateau_);
0502   bool pTOnPlateau = (pt > pTCutForPlateau_);
0503 
0504   if (dxyOnPlateau || !requirePlateau) {
0505     (mes.h_pt)->Fill(pt);
0506   }
0507 
0508   if ((pTOnPlateau && dxyOnPlateau) || !requirePlateau) {
0509     (mes.h_eta)->Fill(eta);
0510     (mes.h_phi)->Fill(phi);
0511     (mes.h_dz)->Fill(dz);
0512     (mes.h_dzWRTpv)->Fill(dzWRTpv);
0513     (mes.h_charge)->Fill(charge);
0514     (mes.h_hits)->Fill(nhits);
0515     (mes.h_onlinelumi)->Fill(onlinelumi);
0516     (mes.h_ls)->Fill(ls);
0517     (mes.h_PU)->Fill(PU);
0518   }
0519 
0520   if (pTOnPlateau || !requirePlateau) {
0521     (mes.h_dxy)->Fill(dxy);
0522     (mes.h_dxyWRTpv)->Fill(dxyWRTpv);
0523   }
0524 
0525   (mes.h_pt_vs_eta)->Fill(eta, pt);
0526 }
0527 
0528 void TrackToTrackComparisonHists::fill_matching_tracks_histos(
0529     matchingME& mes, reco::Track* mon, reco::Track* ref, reco::BeamSpot* bs, reco::Vertex* pv) {
0530   float mon_pt = mon->pt();
0531   float mon_eta = mon->eta();
0532   float mon_phi = mon->phi();
0533   float mon_dxy = mon->dxy(bs->position());
0534   float mon_dz = mon->dz(bs->position());
0535   float mon_dxyWRTpv = mon->dxy(pv->position());
0536   float mon_dzWRTpv = mon->dz(pv->position());
0537   float mon_charge = mon->charge();
0538   float mon_nhits = mon->hitPattern().numberOfValidHits();
0539 
0540   float ref_pt = ref->pt();
0541   float ref_eta = ref->eta();
0542   float ref_phi = ref->phi();
0543   float ref_dxy = ref->dxy(bs->position());
0544   float ref_dz = ref->dz(bs->position());
0545   float ref_dxyWRTpv = ref->dxy(pv->position());
0546   float ref_dzWRTpv = ref->dz(pv->position());
0547   float ref_charge = ref->charge();
0548   float ref_nhits = ref->hitPattern().numberOfValidHits();
0549 
0550   (mes.h_hits_vs_hits)->Fill(ref_nhits, mon_nhits);
0551   (mes.h_pt_vs_pt)->Fill(ref_pt, mon_pt);
0552   (mes.h_eta_vs_eta)->Fill(ref_eta, mon_eta);
0553   (mes.h_phi_vs_phi)->Fill(ref_phi, mon_phi);
0554 
0555   (mes.h_dPt)->Fill(ref_pt - mon_pt);
0556   (mes.h_dEta)->Fill(ref_eta - mon_eta);
0557   (mes.h_dPhi)->Fill(ref_phi - mon_phi);
0558   (mes.h_dDxy)->Fill(ref_dxy - mon_dxy);
0559   (mes.h_dDz)->Fill(ref_dz - mon_dz);
0560   (mes.h_dDxyWRTpv)->Fill(ref_dxyWRTpv - mon_dxyWRTpv);
0561   (mes.h_dDzWRTpv)->Fill(ref_dzWRTpv - mon_dzWRTpv);
0562   (mes.h_dCharge)->Fill(ref_charge - mon_charge);
0563   (mes.h_dHits)->Fill(ref_nhits - mon_nhits);
0564 }
0565 
0566 void TrackToTrackComparisonHists::initialize_parameter(const edm::ParameterSet& iConfig) {
0567   const edm::ParameterSet& pset = iConfig.getParameter<edm::ParameterSet>("histoPSet");
0568 
0569   Eta_rangeMin = pset.getParameter<double>("Eta_rangeMin");
0570   Eta_rangeMax = pset.getParameter<double>("Eta_rangeMax");
0571   Eta_nbin = pset.getParameter<unsigned int>("Eta_nbin");
0572 
0573   Pt_rangeMin = pset.getParameter<double>("Pt_rangeMin");
0574   Pt_rangeMax = pset.getParameter<double>("Pt_rangeMax");
0575   Pt_nbin = pset.getParameter<unsigned int>("Pt_nbin");
0576 
0577   Phi_rangeMin = pset.getParameter<double>("Phi_rangeMin");
0578   Phi_rangeMax = pset.getParameter<double>("Phi_rangeMax");
0579   Phi_nbin = pset.getParameter<unsigned int>("Phi_nbin");
0580 
0581   Dxy_rangeMin = pset.getParameter<double>("Dxy_rangeMin");
0582   Dxy_rangeMax = pset.getParameter<double>("Dxy_rangeMax");
0583   Dxy_nbin = pset.getParameter<unsigned int>("Dxy_nbin");
0584 
0585   Dz_rangeMin = pset.getParameter<double>("Dz_rangeMin");
0586   Dz_rangeMax = pset.getParameter<double>("Dz_rangeMax");
0587   Dz_nbin = pset.getParameter<unsigned int>("Dz_nbin");
0588 
0589   ptRes_rangeMin = pset.getParameter<double>("ptRes_rangeMin");
0590   ptRes_rangeMax = pset.getParameter<double>("ptRes_rangeMax");
0591   ptRes_nbin = pset.getParameter<unsigned int>("ptRes_nbin");
0592 
0593   phiRes_rangeMin = pset.getParameter<double>("phiRes_rangeMin");
0594   phiRes_rangeMax = pset.getParameter<double>("phiRes_rangeMax");
0595   phiRes_nbin = pset.getParameter<unsigned int>("phiRes_nbin");
0596 
0597   etaRes_rangeMin = pset.getParameter<double>("etaRes_rangeMin");
0598   etaRes_rangeMax = pset.getParameter<double>("etaRes_rangeMax");
0599   etaRes_nbin = pset.getParameter<unsigned int>("etaRes_nbin");
0600 
0601   dxyRes_rangeMin = pset.getParameter<double>("dxyRes_rangeMin");
0602   dxyRes_rangeMax = pset.getParameter<double>("dxyRes_rangeMax");
0603   dxyRes_nbin = pset.getParameter<unsigned int>("dxyRes_nbin");
0604 
0605   dzRes_rangeMin = pset.getParameter<double>("dzRes_rangeMin");
0606   dzRes_rangeMax = pset.getParameter<double>("dzRes_rangeMax");
0607   dzRes_nbin = pset.getParameter<unsigned int>("dzRes_nbin");
0608 
0609   ls_rangeMin = pset.getParameter<unsigned int>("ls_rangeMin");
0610   ls_rangeMax = pset.getParameter<unsigned int>("ls_rangeMax");
0611   ls_nbin = pset.getParameter<unsigned int>("ls_nbin");
0612 
0613   onlinelumi_rangeMin = pset.getParameter<double>("onlinelumi_rangeMin");
0614   onlinelumi_rangeMax = pset.getParameter<double>("onlinelumi_rangeMax");
0615   onlinelumi_nbin = pset.getParameter<unsigned int>("onlinelumi_nbin");
0616 
0617   PU_rangeMin = pset.getParameter<double>("PU_rangeMin");
0618   PU_rangeMax = pset.getParameter<double>("PU_rangeMax");
0619   PU_nbin = pset.getParameter<unsigned int>("PU_nbin");
0620 }
0621 
0622 void TrackToTrackComparisonHists::fillHistoPSetDescription(edm::ParameterSetDescription& pset) {
0623   pset.add<double>("Eta_rangeMin", -2.5);
0624   pset.add<double>("Eta_rangeMax", 2.5);
0625   pset.add<unsigned int>("Eta_nbin", 50);
0626 
0627   pset.add<double>("Pt_rangeMin", 0.1);
0628   pset.add<double>("Pt_rangeMax", 100.0);
0629   pset.add<unsigned int>("Pt_nbin", 1000);
0630 
0631   pset.add<double>("Phi_rangeMin", -3.1416);
0632   pset.add<double>("Phi_rangeMax", 3.1416);
0633   pset.add<unsigned int>("Phi_nbin", 36);
0634 
0635   pset.add<double>("Dxy_rangeMin", -1.0);
0636   pset.add<double>("Dxy_rangeMax", 1.0);
0637   pset.add<unsigned int>("Dxy_nbin", 300);
0638 
0639   pset.add<double>("Dz_rangeMin", -30.0);
0640   pset.add<double>("Dz_rangeMax", 30.0);
0641   pset.add<unsigned int>("Dz_nbin", 60);
0642 
0643   pset.add<double>("ptRes_rangeMin", -0.1);
0644   pset.add<double>("ptRes_rangeMax", 0.1);
0645   pset.add<unsigned int>("ptRes_nbin", 100);
0646 
0647   pset.add<double>("phiRes_rangeMin", -0.01);
0648   pset.add<double>("phiRes_rangeMax", 0.01);
0649   pset.add<unsigned int>("phiRes_nbin", 300);
0650 
0651   pset.add<double>("etaRes_rangeMin", -0.01);
0652   pset.add<double>("etaRes_rangeMax", 0.01);
0653   pset.add<unsigned int>("etaRes_nbin", 300);
0654 
0655   pset.add<double>("dxyRes_rangeMin", -0.05);
0656   pset.add<double>("dxyRes_rangeMax", 0.05);
0657   pset.add<unsigned int>("dxyRes_nbin", 500);
0658 
0659   pset.add<double>("dzRes_rangeMin", -0.05);
0660   pset.add<double>("dzRes_rangeMax", 0.05);
0661   pset.add<unsigned int>("dzRes_nbin", 150);
0662 
0663   pset.add<unsigned int>("ls_rangeMin", 0);
0664   pset.add<unsigned int>("ls_rangeMax", 3000);
0665   pset.add<unsigned int>("ls_nbin", 300);
0666 
0667   pset.add<double>("onlinelumi_rangeMin", 0.0);
0668   pset.add<double>("onlinelumi_rangeMax", 20000.0);
0669   pset.add<unsigned int>("onlinelumi_nbin", 200);
0670 
0671   pset.add<double>("PU_rangeMin", 0.0);
0672   pset.add<double>("PU_rangeMax", 120.0);
0673   pset.add<unsigned int>("PU_nbin", 120);
0674 }