Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-04-08 03:34:24

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