Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-04-17 02:46:42

0001 /*
0002  *  See header file for a description of this class.
0003  **
0004  *  \author Suchandra Dutta , Giorgia Mila
0005  */
0006 
0007 #include "DataFormats/TrackReco/interface/Track.h"
0008 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0009 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0010 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0011 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0012 
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include "FWCore/Utilities/interface/InputTag.h"
0015 #include "FWCore/Utilities/interface/transform.h"
0016 #include "DQMServices/Core/interface/DQMStore.h"
0017 #include "DQM/TrackingMonitor/interface/TrackBuildingAnalyzer.h"
0018 #include "DQM/TrackingMonitor/interface/TrackAnalyzer.h"
0019 #include "DQM/TrackingMonitor/interface/VertexMonitor.h"
0020 #include "DQM/TrackingMonitor/interface/TrackingMonitor.h"
0021 
0022 #include "FWCore/ParameterSet/interface/Registry.h"
0023 
0024 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0025 #include "TrackingTools/PatternTools/interface/TSCBLBuilderNoMaterial.h"
0026 #include "TrackingTools/PatternTools/interface/TSCPBuilderNoMaterial.h"
0027 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0028 #include "CommonTools/TriggerUtils/interface/GenericTriggerEventFlag.h"
0029 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0030 
0031 #include "DataFormats/Candidate/interface/Candidate.h"
0032 #include "DataFormats/Candidate/interface/CandidateFwd.h"
0033 
0034 #include "DQM/TrackingMonitor/interface/GetLumi.h"
0035 
0036 #include <string>
0037 
0038 #ifdef VI_DEBUG
0039 #define COUT(x) std::cout << x << ' '
0040 #else
0041 #define COUT(x) LogDebug(x)
0042 #endif
0043 
0044 // TrackingMonitor
0045 // ----------------------------------------------------------------------------------//
0046 
0047 TrackingMonitor::TrackingMonitor(const edm::ParameterSet& iConfig)
0048     : confID_(iConfig.id()),
0049       theTrackBuildingAnalyzer(new TrackBuildingAnalyzer(iConfig)),
0050       NumberOfTracks(nullptr),
0051       NumberOfMeanRecHitsPerTrack(nullptr),
0052       NumberOfMeanLayersPerTrack(nullptr)
0053       //    , NumberOfGoodTracks(NULL)
0054       ,
0055       FractionOfGoodTracks(nullptr),
0056       NumberOfTrackingRegions(nullptr),
0057       NumberOfSeeds(nullptr),
0058       NumberOfSeeds_lumiFlag(nullptr),
0059       NumberOfTrackCandidates(nullptr),
0060       FractionCandidatesOverSeeds(nullptr)
0061       //    , NumberOfGoodTrkVsClus(NULL)
0062       ,
0063       NumberEventsOfVsLS(nullptr),
0064       NumberOfTracksVsLS(nullptr)
0065       //    , NumberOfGoodTracksVsLS(NULL)
0066       ,
0067       GoodTracksFractionVsLS(nullptr)
0068       //    , GoodTracksNumberOfRecHitsPerTrackVsLS(NULL)
0069       // ADD by Mia for PU monitoring
0070       // vertex plots to be moved in ad hoc class
0071       ,
0072       NumberOfGoodPVtxVsLS(nullptr),
0073       NumberOfGoodPVtxWO0VsLS(nullptr),
0074       NumberEventsOfVsBX(nullptr),
0075       NumberOfTracksVsBX(nullptr),
0076       GoodTracksFractionVsBX(nullptr),
0077       NumberOfRecHitsPerTrackVsBX(nullptr),
0078       NumberOfGoodPVtxVsBX(nullptr),
0079       NumberOfGoodPVtxWO0VsBX(nullptr),
0080       NumberOfTracksVsBXlumi(nullptr),
0081       NumberOfTracksVsGoodPVtx(nullptr),
0082       NumberOfTracksVsPUPVtx(nullptr),
0083       NumberEventsOfVsGoodPVtx(nullptr),
0084       GoodTracksFractionVsGoodPVtx(nullptr),
0085       NumberOfRecHitsPerTrackVsGoodPVtx(nullptr),
0086       NumberOfPVtxVsGoodPVtx(nullptr),
0087       NumberOfPixelClustersVsGoodPVtx(nullptr),
0088       NumberOfStripClustersVsGoodPVtx(nullptr),
0089       NumberEventsOfVsLUMI(nullptr),
0090       NumberOfTracksVsLUMI(nullptr),
0091       GoodTracksFractionVsLUMI(nullptr),
0092       NumberOfRecHitsPerTrackVsLUMI(nullptr),
0093       NumberOfGoodPVtxVsLUMI(nullptr),
0094       NumberOfGoodPVtxWO0VsLUMI(nullptr),
0095       NumberOfPixelClustersVsLUMI(nullptr),
0096       NumberOfStripClustersVsLUMI(nullptr),
0097       NumberOfTracks_lumiFlag(nullptr)
0098       //    , NumberOfGoodTracks_lumiFlag(NULL)
0099 
0100       ,
0101       builderName(iConfig.getParameter<std::string>("TTRHBuilder")),
0102       doTrackerSpecific_(iConfig.getParameter<bool>("doTrackerSpecific")),
0103       doLumiAnalysis(iConfig.getParameter<bool>("doLumiAnalysis")),
0104       doProfilesVsLS_(iConfig.getParameter<bool>("doProfilesVsLS")),
0105       doAllSeedPlots(iConfig.getParameter<bool>("doSeedParameterHistos")),
0106       doAllPlots(iConfig.getParameter<bool>("doAllPlots")),
0107       doGeneralPropertiesPlots_(iConfig.getParameter<bool>("doGeneralPropertiesPlots")),
0108       doHitPropertiesPlots_(iConfig.getParameter<bool>("doHitPropertiesPlots")),
0109       doTkCandPlots(iConfig.getParameter<bool>("doTrackCandHistos")),
0110       doPUmonitoring_(iConfig.getParameter<bool>("doPUmonitoring")),
0111       genTriggerEventFlag_(new GenericTriggerEventFlag(
0112           iConfig.getParameter<edm::ParameterSet>("genericTriggerEventPSet"), consumesCollector(), *this)),
0113       numSelection_(iConfig.getParameter<std::string>("numCut")),
0114       denSelection_(iConfig.getParameter<std::string>("denCut")),
0115       pvNDOF_(iConfig.getParameter<int>("pvNDOF")),
0116       forceSCAL_(iConfig.getParameter<bool>("forceSCAL")) {
0117   edm::ConsumesCollector c{consumesCollector()};
0118   theTrackAnalyzer = new tadqm::TrackAnalyzer(iConfig, c);
0119 
0120   // input tags for collections from the configuration
0121   bsSrc_ = iConfig.getParameter<edm::InputTag>("beamSpot");
0122   pvSrc_ = iConfig.getParameter<edm::InputTag>("primaryVertex");
0123   bsSrcToken_ = consumes<reco::BeamSpot>(bsSrc_);
0124   pvSrcToken_ = mayConsume<reco::VertexCollection>(pvSrc_);
0125 
0126   lumiscalersToken_ = consumes<LumiScalersCollection>(iConfig.getParameter<edm::InputTag>("scal"));
0127   metaDataToken_ = consumes<OnlineLuminosityRecord>(iConfig.getParameter<edm::InputTag>("metadata"));
0128 
0129   edm::InputTag alltrackProducer = iConfig.getParameter<edm::InputTag>("allTrackProducer");
0130   edm::InputTag trackProducer = iConfig.getParameter<edm::InputTag>("TrackProducer");
0131   edm::InputTag tcProducer = iConfig.getParameter<edm::InputTag>("TCProducer");
0132   edm::InputTag seedProducer = iConfig.getParameter<edm::InputTag>("SeedProducer");
0133   allTrackToken_ = consumes<edm::View<reco::Track> >(alltrackProducer);
0134   trackToken_ = consumes<edm::View<reco::Track> >(trackProducer);
0135   trackCandidateToken_ = consumes<TrackCandidateCollection>(tcProducer);
0136   seedToken_ = consumes<edm::View<TrajectorySeed> >(seedProducer);
0137   seedStopInfoToken_ = consumes<std::vector<SeedStopInfo> >(tcProducer);
0138 
0139   doMVAPlots = iConfig.getParameter<bool>("doMVAPlots");
0140   if (doMVAPlots) {
0141     mvaQualityTokens_ = edm::vector_transform(
0142         iConfig.getParameter<std::vector<std::string> >("MVAProducers"), [&](const std::string& tag) {
0143           return std::make_tuple(consumes<MVACollection>(edm::InputTag(tag, "MVAValues")),
0144                                  consumes<QualityMaskCollection>(edm::InputTag(tag, "QualityMasks")));
0145         });
0146     mvaTrackToken_ = consumes<edm::View<reco::Track> >(iConfig.getParameter<edm::InputTag>("TrackProducerForMVA"));
0147   }
0148 
0149   doRegionPlots = iConfig.getParameter<bool>("doRegionPlots");
0150   doRegionCandidatePlots = iConfig.getParameter<bool>("doRegionCandidatePlots");
0151   if (doRegionPlots) {
0152     const auto& regionTag = iConfig.getParameter<edm::InputTag>("RegionProducer");
0153     if (!regionTag.label().empty()) {
0154       regionToken_ = consumes<edm::OwnVector<TrackingRegion> >(regionTag);
0155     }
0156     const auto& regionLayersTag = iConfig.getParameter<edm::InputTag>("RegionSeedingLayersProducer");
0157     if (!regionLayersTag.label().empty()) {
0158       if (!regionToken_.isUninitialized()) {
0159         throw cms::Exception("Configuration")
0160             << "Only one of 'RegionProducer' and 'RegionSeedingLayersProducer' can be non-empty, now both are.";
0161       }
0162       regionLayerSetsToken_ = consumes<TrackingRegionsSeedingLayerSets>(regionLayersTag);
0163     } else if (regionToken_.isUninitialized()) {
0164       throw cms::Exception("Configuration") << "With doRegionPlots=True either 'RegionProducer' or "
0165                                                "'RegionSeedingLayersProducer' must be non-empty, now both are empty.";
0166     }
0167 
0168     if (doRegionCandidatePlots) {
0169       regionCandidateToken_ = consumes<reco::CandidateView>(iConfig.getParameter<edm::InputTag>("RegionCandidates"));
0170     }
0171   }
0172 
0173   edm::InputTag stripClusterInputTag_ = iConfig.getParameter<edm::InputTag>("stripCluster");
0174   edm::InputTag pixelClusterInputTag_ = iConfig.getParameter<edm::InputTag>("pixelCluster");
0175   stripClustersToken_ = mayConsume<edmNew::DetSetVector<SiStripCluster> >(stripClusterInputTag_);
0176   pixelClustersToken_ = mayConsume<edmNew::DetSetVector<SiPixelCluster> >(pixelClusterInputTag_);
0177 
0178   runTrackBuildingAnalyzerForSeed =
0179       (doAllSeedPlots || iConfig.getParameter<bool>("doSeedPTHisto") || iConfig.getParameter<bool>("doSeedETAHisto") ||
0180        iConfig.getParameter<bool>("doSeedPHIHisto") || iConfig.getParameter<bool>("doSeedPHIVsETAHisto") ||
0181        iConfig.getParameter<bool>("doSeedThetaHisto") || iConfig.getParameter<bool>("doSeedQHisto") ||
0182        iConfig.getParameter<bool>("doSeedDxyHisto") || iConfig.getParameter<bool>("doSeedDzHisto") ||
0183        iConfig.getParameter<bool>("doSeedNRecHitsHisto") || iConfig.getParameter<bool>("doSeedNVsPhiProf") ||
0184        iConfig.getParameter<bool>("doSeedNVsEtaProf"));
0185 
0186   if (doTkCandPlots || doAllSeedPlots || runTrackBuildingAnalyzerForSeed) {
0187     magneticFieldToken_ = esConsumes<MagneticField, IdealMagneticFieldRecord>();
0188     transientTrackingRecHitBuilderToken_ =
0189         esConsumes<TransientTrackingRecHitBuilder, TransientRecHitRecord>(edm::ESInputTag("", builderName));
0190   }
0191 
0192   doFractionPlot_ = true;
0193   if (alltrackProducer.label() == trackProducer.label())
0194     doFractionPlot_ = false;
0195 
0196   Quality_ = iConfig.getParameter<std::string>("Quality");
0197   AlgoName_ = iConfig.getParameter<std::string>("AlgoName");
0198 
0199   // get flag from the configuration
0200   doPlotsVsBXlumi_ = iConfig.getParameter<bool>("doPlotsVsBXlumi");
0201   if (doPlotsVsBXlumi_)
0202     theLumiDetails_ = new GetLumi(iConfig.getParameter<edm::ParameterSet>("BXlumiSetup"), c);
0203   doPlotsVsGoodPVtx_ = iConfig.getParameter<bool>("doPlotsVsGoodPVtx");
0204   doPlotsVsLUMI_ = iConfig.getParameter<bool>("doPlotsVsLUMI");
0205   doPlotsVsBX_ = iConfig.getParameter<bool>("doPlotsVsBX");
0206 
0207   if (doPUmonitoring_) {
0208     std::vector<edm::InputTag> primaryVertexInputTags =
0209         iConfig.getParameter<std::vector<edm::InputTag> >("primaryVertexInputTags");
0210     std::vector<edm::InputTag> selPrimaryVertexInputTags =
0211         iConfig.getParameter<std::vector<edm::InputTag> >("selPrimaryVertexInputTags");
0212     std::vector<std::string> pvLabels = iConfig.getParameter<std::vector<std::string> >("pvLabels");
0213 
0214     if (primaryVertexInputTags.size() == pvLabels.size() and
0215         primaryVertexInputTags.size() == selPrimaryVertexInputTags.size()) {
0216       for (size_t i = 0; i < primaryVertexInputTags.size(); i++) {
0217         edm::InputTag iPVinputTag = primaryVertexInputTags[i];
0218         edm::InputTag iSelPVinputTag = selPrimaryVertexInputTags[i];
0219         std::string iPVlabel = pvLabels[i];
0220 
0221         theVertexMonitor.push_back(new VertexMonitor(iConfig, iPVinputTag, iSelPVinputTag, iPVlabel, c));
0222       }
0223     }
0224   }
0225 }
0226 
0227 TrackingMonitor::~TrackingMonitor() {
0228   if (theTrackAnalyzer)
0229     delete theTrackAnalyzer;
0230   if (theTrackBuildingAnalyzer)
0231     delete theTrackBuildingAnalyzer;
0232   if (doPUmonitoring_)
0233     for (size_t i = 0; i < theVertexMonitor.size(); i++)
0234       if (theVertexMonitor[i])
0235         delete theVertexMonitor[i];
0236   if (genTriggerEventFlag_)
0237     delete genTriggerEventFlag_;
0238 }
0239 
0240 void TrackingMonitor::bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& iRun, edm::EventSetup const& iSetup) {
0241   // parameters from the configuration
0242   auto const* conf = edm::pset::Registry::instance()->getMapped(confID_);
0243   assert(conf != nullptr);
0244   std::string Quality = conf->getParameter<std::string>("Quality");
0245   std::string AlgoName = conf->getParameter<std::string>("AlgoName");
0246   MEFolderName = conf->getParameter<std::string>("FolderName");
0247   std::string Folder = MEFolderName.substr(0, 2);
0248 
0249   // test for the Quality veriable validity
0250   if (!Quality_.empty()) {
0251     if (Quality_ != "highPurity" && Quality_ != "tight" && Quality_ != "loose") {
0252       edm::LogWarning("TrackingMonitor") << "Qualty Name is invalid, using no quality criterea by default";
0253       Quality_ = "";
0254     }
0255   }
0256 
0257   // use the AlgoName and Quality Name
0258   std::string CategoryName = !Quality_.empty() ? AlgoName_ + "_" + Quality_ : AlgoName_;
0259 
0260   // get binning from the configuration
0261   int TKNoBin = conf->getParameter<int>("TkSizeBin");
0262   double TKNoMin = conf->getParameter<double>("TkSizeMin");
0263   double TKNoMax = conf->getParameter<double>("TkSizeMax");
0264 
0265   int TCNoBin = conf->getParameter<int>("TCSizeBin");
0266   double TCNoMin = conf->getParameter<double>("TCSizeMin");
0267   double TCNoMax = conf->getParameter<double>("TCSizeMax");
0268 
0269   int TKNoSeedBin = conf->getParameter<int>("TkSeedSizeBin");
0270   double TKNoSeedMin = conf->getParameter<double>("TkSeedSizeMin");
0271   double TKNoSeedMax = conf->getParameter<double>("TkSeedSizeMax");
0272 
0273   //int RecHitBin = conf->getParameter<int>("RecHitBin");
0274   double RecHitMin = conf->getParameter<double>("RecHitMin");
0275   double RecHitMax = conf->getParameter<double>("RecHitMax");
0276 
0277   int MeanHitBin = conf->getParameter<int>("MeanHitBin");
0278   double MeanHitMin = conf->getParameter<double>("MeanHitMin");
0279   double MeanHitMax = conf->getParameter<double>("MeanHitMax");
0280 
0281   int MeanLayBin = conf->getParameter<int>("MeanLayBin");
0282   double MeanLayMin = conf->getParameter<double>("MeanLayMin");
0283   double MeanLayMax = conf->getParameter<double>("MeanLayMax");
0284 
0285   int LSBin = conf->getParameter<int>("LSBin");
0286   int LSMin = conf->getParameter<double>("LSMin");
0287   int LSMax = conf->getParameter<double>("LSMax");
0288 
0289   std::string StateName = conf->getParameter<std::string>("MeasurementState");
0290   if (StateName != "OuterSurface" && StateName != "InnerSurface" && StateName != "ImpactPoint" &&
0291       StateName != "default" && StateName != "All") {
0292     // print warning
0293     edm::LogWarning("TrackingMonitor") << "State Name is invalid, using 'ImpactPoint' by default";
0294   }
0295 
0296   ibooker.setCurrentFolder(MEFolderName);
0297 
0298   // book the General Property histograms
0299   // ---------------------------------------------------------------------------------//
0300 
0301   if (doGeneralPropertiesPlots_ || doAllPlots) {
0302     ibooker.setCurrentFolder(MEFolderName + "/GeneralProperties");
0303 
0304     histname = "NumberOfTracks_" + CategoryName;
0305     // MODIFY by Mia in order to cope w/ high multiplicity
0306     NumberOfTracks = ibooker.book1D(histname, histname, TKNoBin, TKNoMin, TKNoMax);
0307     NumberOfTracks->setAxisTitle("Number of Tracks per Event", 1);
0308     NumberOfTracks->setAxisTitle("Number of Events", 2);
0309 
0310     if (Folder == "Tr") {
0311       histname = "NumberOfTracks_PUvtx_" + CategoryName;
0312       NumberOfTracks_PUvtx = ibooker.book1D(histname, histname, TKNoBin, TKNoMin, TKNoMax);
0313       NumberOfTracks_PUvtx->setAxisTitle("Number of Tracks per Event (matched a PU vertex)", 1);
0314       NumberOfTracks_PUvtx->setAxisTitle("Number of Events", 2);
0315 
0316       histname = "NumberofTracks_Hardvtx_" + CategoryName;
0317       NumberofTracks_Hardvtx =
0318           ibooker.book1D(histname, histname, TKNoBin / 10, TKNoMin, ((TKNoMax - TKNoMin) / 10) - 0.5);
0319       NumberofTracks_Hardvtx->setAxisTitle("Number of Tracks per Event (matched main vertex)", 1);
0320       NumberofTracks_Hardvtx->setAxisTitle("Number of Events", 2);
0321 
0322       histname = "NumberofTracks_Hardvtx_PUvtx_" + CategoryName;
0323       NumberofTracks_Hardvtx_PUvtx = ibooker.book1D(histname, histname, 2, 0., 2.);
0324       NumberofTracks_Hardvtx_PUvtx->setAxisTitle("Number of Tracks per PU/Hard vertex", 1);
0325       NumberofTracks_Hardvtx_PUvtx->setAxisTitle("Number of Tracks", 2);
0326       NumberofTracks_Hardvtx_PUvtx->setBinLabel(1, "PU_Vertex");
0327       NumberofTracks_Hardvtx_PUvtx->setBinLabel(2, "Hard_Vertex");
0328     }
0329 
0330     histname = "NumberOfMeanRecHitsPerTrack_" + CategoryName;
0331     NumberOfMeanRecHitsPerTrack = ibooker.book1D(histname, histname, MeanHitBin, MeanHitMin, MeanHitMax);
0332     NumberOfMeanRecHitsPerTrack->setAxisTitle("Mean number of valid RecHits per Track", 1);
0333     NumberOfMeanRecHitsPerTrack->setAxisTitle("Entries", 2);
0334 
0335     histname = "NumberOfMeanLayersPerTrack_" + CategoryName;
0336     NumberOfMeanLayersPerTrack = ibooker.book1D(histname, histname, MeanLayBin, MeanLayMin, MeanLayMax);
0337     NumberOfMeanLayersPerTrack->setAxisTitle("Mean number of Layers per Track", 1);
0338     NumberOfMeanLayersPerTrack->setAxisTitle("Entries", 2);
0339 
0340     if (doFractionPlot_) {
0341       histname = "FractionOfGoodTracks_" + CategoryName;
0342       FractionOfGoodTracks = ibooker.book1D(histname, histname, 101, -0.005, 1.005);
0343       FractionOfGoodTracks->setAxisTitle("Fraction of Tracks (w.r.t. generalTracks)", 1);
0344       FractionOfGoodTracks->setAxisTitle("Entries", 2);
0345     }
0346   }
0347 
0348   if (doLumiAnalysis) {
0349     // add by Mia in order to deal with LS transitions
0350     ibooker.setCurrentFolder(MEFolderName + "/LSanalysis");
0351     auto scope = DQMStore::IBooker::UseLumiScope(ibooker);
0352 
0353     histname = "NumberOfTracks_lumiFlag_" + CategoryName;
0354     NumberOfTracks_lumiFlag = ibooker.book1D(histname, histname, TKNoBin, TKNoMin, TKNoMax);
0355     NumberOfTracks_lumiFlag->setAxisTitle("Number of Tracks per Event", 1);
0356     NumberOfTracks_lumiFlag->setAxisTitle("Number of Events", 2);
0357   }
0358 
0359   // book profile plots vs LS :
0360   //---------------------------
0361 
0362   if (doProfilesVsLS_ || doAllPlots) {
0363     ibooker.setCurrentFolder(MEFolderName + "/GeneralProperties");
0364 
0365     histname = "NumberOfTracksVsLS_" + CategoryName;
0366     NumberOfTracksVsLS = ibooker.bookProfile(histname, histname, LSBin, LSMin, LSMax, TKNoMin, TKNoMax, "");
0367     NumberOfTracksVsLS->getTH1()->SetCanExtend(TH1::kAllAxes);
0368     NumberOfTracksVsLS->setAxisTitle("#Lumi section", 1);
0369     NumberOfTracksVsLS->setAxisTitle("Number of  Tracks", 2);
0370 
0371     histname = "NumberOfRecHitsPerTrackVsLS_" + CategoryName;
0372     NumberOfRecHitsPerTrackVsLS =
0373         ibooker.bookProfile(histname, histname, LSBin, LSMin, LSMax, RecHitMin, RecHitMax * 5, "");
0374     NumberOfRecHitsPerTrackVsLS->getTH1()->SetCanExtend(TH1::kAllAxes);
0375     NumberOfRecHitsPerTrackVsLS->setAxisTitle("#Lumi section", 1);
0376     NumberOfRecHitsPerTrackVsLS->setAxisTitle("Mean number of Valid RecHits per track", 2);
0377 
0378     histname = "NumberEventsVsLS_" + CategoryName;
0379     NumberEventsOfVsLS = ibooker.book1D(histname, histname, LSBin, LSMin, LSMax);
0380     NumberEventsOfVsLS->getTH1()->SetCanExtend(TH1::kAllAxes);
0381     NumberEventsOfVsLS->setAxisTitle("#Lumi section", 1);
0382     NumberEventsOfVsLS->setAxisTitle("Number of events", 2);
0383 
0384     edm::ParameterSet ParametersGoodPVtx = conf->getParameter<edm::ParameterSet>("GoodPVtx");
0385 
0386     double GoodPVtxMin = ParametersGoodPVtx.getParameter<double>("GoodPVtxMin");
0387     double GoodPVtxMax = ParametersGoodPVtx.getParameter<double>("GoodPVtxMax");
0388 
0389     histname = "NumberOfGoodPVtxVsLS_" + CategoryName;
0390     NumberOfGoodPVtxVsLS =
0391         ibooker.bookProfile(histname, histname, LSBin, LSMin, LSMax, GoodPVtxMin, 3. * GoodPVtxMax, "");
0392     NumberOfGoodPVtxVsLS->getTH1()->SetCanExtend(TH1::kAllAxes);
0393     NumberOfGoodPVtxVsLS->setAxisTitle("#Lumi section", 1);
0394     NumberOfGoodPVtxVsLS->setAxisTitle("Mean number of good PV", 2);
0395 
0396     histname = "NumberOfGoodPVtxWO0VsLS_" + CategoryName;
0397     NumberOfGoodPVtxWO0VsLS =
0398         ibooker.bookProfile(histname, histname, LSBin, LSMin, LSMax, GoodPVtxMin, 3. * GoodPVtxMax, "");
0399     NumberOfGoodPVtxWO0VsLS->setAxisTitle("#Lumi section", 1);
0400     NumberOfGoodPVtxWO0VsLS->setAxisTitle("Mean number of good PV", 2);
0401 
0402     if (doFractionPlot_) {
0403       histname = "GoodTracksFractionVsLS_" + CategoryName;
0404       GoodTracksFractionVsLS = ibooker.bookProfile(histname, histname, LSBin, LSMin, LSMax, 0, 1.1, "");
0405       GoodTracksFractionVsLS->getTH1()->SetCanExtend(TH1::kAllAxes);
0406       GoodTracksFractionVsLS->setAxisTitle("#Lumi section", 1);
0407       GoodTracksFractionVsLS->setAxisTitle("Fraction of Good Tracks", 2);
0408     }
0409 
0410     if (doPlotsVsBX_ || doAllPlots) {
0411       ibooker.setCurrentFolder(MEFolderName + "/BXanalysis");
0412       int BXBin = 3564;
0413       double BXMin = 0.5;
0414       double BXMax = 3564.5;
0415 
0416       histname = "NumberEventsVsBX_" + CategoryName;
0417       NumberEventsOfVsBX = ibooker.book1D(histname, histname, BXBin, BXMin, BXMax);
0418       NumberEventsOfVsBX->setAxisTitle("BX", 1);
0419       NumberEventsOfVsBX->setAxisTitle("Number of events", 2);
0420 
0421       histname = "NumberOfTracksVsBX_" + CategoryName;
0422       NumberOfTracksVsBX = ibooker.bookProfile(histname, histname, BXBin, BXMin, BXMax, TKNoMin, TKNoMax * 3., "");
0423       NumberOfTracksVsBX->setAxisTitle("BX", 1);
0424       NumberOfTracksVsBX->setAxisTitle("Number of  Tracks", 2);
0425 
0426       histname = "NumberOfRecHitsPerTrackVsBX_" + CategoryName;
0427       NumberOfRecHitsPerTrackVsBX =
0428           ibooker.bookProfile(histname, histname, BXBin, BXMin, BXMax, RecHitMin, RecHitMax * 5, "");
0429       NumberOfRecHitsPerTrackVsBX->setAxisTitle("BX", 1);
0430       NumberOfRecHitsPerTrackVsBX->setAxisTitle("Mean number of Valid RecHits per track", 2);
0431 
0432       histname = "NumberOfGoodPVtxVsBX_" + CategoryName;
0433       NumberOfGoodPVtxVsBX =
0434           ibooker.bookProfile(histname, histname, BXBin, BXMin, BXMax, GoodPVtxMin, 3. * GoodPVtxMax, "");
0435       NumberOfGoodPVtxVsBX->setAxisTitle("BX", 1);
0436       NumberOfGoodPVtxVsBX->setAxisTitle("Mean number of good PV", 2);
0437 
0438       histname = "NumberOfGoodPVtxWO0VsBX_" + CategoryName;
0439       NumberOfGoodPVtxWO0VsBX =
0440           ibooker.bookProfile(histname, histname, BXBin, BXMin, BXMax, GoodPVtxMin, 3. * GoodPVtxMax, "");
0441       NumberOfGoodPVtxWO0VsBX->setAxisTitle("BX", 1);
0442       NumberOfGoodPVtxWO0VsBX->setAxisTitle("Mean number of good PV", 2);
0443 
0444       if (doFractionPlot_) {
0445         histname = "GoodTracksFractionVsBX_" + CategoryName;
0446         GoodTracksFractionVsBX = ibooker.bookProfile(histname, histname, BXBin, BXMin, BXMax, 0, 1.1, "");
0447         GoodTracksFractionVsBX->setAxisTitle("BX", 1);
0448         GoodTracksFractionVsBX->setAxisTitle("Fraction of Good Tracks", 2);
0449       }
0450     }
0451   }
0452 
0453   // book PU monitoring plots :
0454   //---------------------------
0455 
0456   if (doPUmonitoring_) {
0457     for (size_t i = 0; i < theVertexMonitor.size(); i++)
0458       theVertexMonitor[i]->initHisto(ibooker);
0459   }
0460 
0461   if (doPlotsVsGoodPVtx_) {
0462     ibooker.setCurrentFolder(MEFolderName + "/PUmonitoring");
0463     // get binning from the configuration
0464     int PVBin = conf->getParameter<int>("PVBin");
0465     float PVMin = conf->getParameter<double>("PVMin");
0466     float PVMax = conf->getParameter<double>("PVMax");
0467 
0468     histname = "NumberOfTracksVsGoodPVtx";
0469     NumberOfTracksVsGoodPVtx = ibooker.bookProfile(histname, histname, PVBin, PVMin, PVMax, TKNoMin, TKNoMax * 5., "");
0470     NumberOfTracksVsGoodPVtx->setAxisTitle("Number of PV", 1);
0471     NumberOfTracksVsGoodPVtx->setAxisTitle("Mean number of Tracks per Event", 2);
0472 
0473     histname = "NumberOfTracksVsPUPVtx";
0474     NumberOfTracksVsPUPVtx = ibooker.bookProfile(histname, histname, PVBin, PVMin, PVMax, 0., TKNoMax * 5., "");
0475     NumberOfTracksVsPUPVtx->setAxisTitle("Number of PU", 1);
0476     NumberOfTracksVsPUPVtx->setAxisTitle("Mean number of Tracks per PUvtx", 2);
0477 
0478     histname = "NumberEventsVsGoodPVtx";
0479     NumberEventsOfVsGoodPVtx = ibooker.book1D(histname, histname, PVBin, PVMin, PVMax);
0480     NumberEventsOfVsGoodPVtx->setAxisTitle("Number of good PV (PU)", 1);
0481     NumberEventsOfVsGoodPVtx->setAxisTitle("Number of events", 2);
0482 
0483     if (doFractionPlot_) {
0484       histname = "GoodTracksFractionVsGoodPVtx";
0485       GoodTracksFractionVsGoodPVtx = ibooker.bookProfile(histname, histname, PVBin, PVMin, PVMax, 0., 1.1, "");
0486       GoodTracksFractionVsGoodPVtx->setAxisTitle("Number of good PV (PU)", 1);
0487       GoodTracksFractionVsGoodPVtx->setAxisTitle("Mean fraction of good tracks", 2);
0488     }
0489 
0490     histname = "NumberOfRecHitsPerTrackVsGoodPVtx";
0491     NumberOfRecHitsPerTrackVsGoodPVtx = ibooker.bookProfile(histname, histname, PVBin, PVMin, PVMax, 0., 200., "");
0492     NumberOfRecHitsPerTrackVsGoodPVtx->setAxisTitle("Number of good PV (PU)", 1);
0493     NumberOfRecHitsPerTrackVsGoodPVtx->setAxisTitle("Mean number of valid rechits per Tracks", 2);
0494 
0495     histname = "NumberOfPVtxVsGoodPVtx";
0496     NumberOfPVtxVsGoodPVtx = ibooker.bookProfile(histname, histname, PVBin, PVMin, PVMax, 0., 3. * PVMax, "");
0497     NumberOfPVtxVsGoodPVtx->setAxisTitle("Number of good PV (PU)", 1);
0498     NumberOfPVtxVsGoodPVtx->setAxisTitle("Mean number of vertices", 2);
0499 
0500     double NClusPxMin = conf->getParameter<double>("NClusPxMin");
0501     double NClusPxMax = conf->getParameter<double>("NClusPxMax");
0502     histname = "NumberOfPixelClustersVsGoodPVtx";
0503     NumberOfPixelClustersVsGoodPVtx =
0504         ibooker.bookProfile(histname, histname, PVBin, PVMin, PVMax, NClusPxMin, 3. * NClusPxMax, "");
0505     NumberOfPixelClustersVsGoodPVtx->setAxisTitle("Number of good PV (PU)", 1);
0506     NumberOfPixelClustersVsGoodPVtx->setAxisTitle("Mean number of pixel clusters", 2);
0507 
0508     double NClusStrMin = conf->getParameter<double>("NClusStrMin");
0509     double NClusStrMax = conf->getParameter<double>("NClusStrMax");
0510     histname = "NumberOfStripClustersVsGoodPVtx";
0511     NumberOfStripClustersVsGoodPVtx =
0512         ibooker.bookProfile(histname, histname, PVBin, PVMin, PVMax, NClusStrMin, 3. * NClusStrMax, "");
0513     NumberOfStripClustersVsGoodPVtx->setAxisTitle("Number of good PV (PU)", 1);
0514     NumberOfStripClustersVsGoodPVtx->setAxisTitle("Mean number of strip clusters", 2);
0515   }
0516 
0517   if (doPlotsVsLUMI_ || doAllPlots) {
0518     ibooker.setCurrentFolder(MEFolderName + "/LUMIanalysis");
0519     int LUMIBin = conf->getParameter<int>("LUMIBin");
0520     float LUMIMin = conf->getParameter<double>("LUMIMin");
0521     float LUMIMax = conf->getParameter<double>("LUMIMax");
0522 
0523     histname = "NumberEventsVsLUMI";
0524     NumberEventsOfVsLUMI = ibooker.book1D(histname, histname, LUMIBin, LUMIMin, LUMIMax);
0525     NumberEventsOfVsLUMI->setAxisTitle("scal lumi [10e30 Hz cm^{-2}]", 1);
0526     NumberEventsOfVsLUMI->setAxisTitle("Number of events", 2);
0527 
0528     histname = "NumberOfTracksVsLUMI";
0529     NumberOfTracksVsLUMI = ibooker.bookProfile(histname, histname, LUMIBin, LUMIMin, LUMIMax, TKNoMin, TKNoMax * 3, "");
0530     NumberOfTracksVsLUMI->setAxisTitle("scal lumi [10e30 Hz cm^{-2}]", 1);
0531     NumberOfTracksVsLUMI->setAxisTitle("Mean number of vertices", 2);
0532 
0533     if (doFractionPlot_) {
0534       histname = "GoodTracksFractionVsLUMI";
0535       GoodTracksFractionVsLUMI = ibooker.bookProfile(histname, histname, LUMIBin, LUMIMin, LUMIMax, 0., 1.1, "");
0536       GoodTracksFractionVsLUMI->setAxisTitle("scal lumi [10e30 Hz cm^{-2}]", 1);
0537       GoodTracksFractionVsLUMI->setAxisTitle("Mean number of vertices", 2);
0538     }
0539 
0540     histname = "NumberOfRecHitsPerTrackVsLUMI";
0541     NumberOfRecHitsPerTrackVsLUMI =
0542         ibooker.bookProfile(histname, histname, LUMIBin, LUMIMin, LUMIMax, RecHitMin, RecHitMax * 5, "");
0543     NumberOfRecHitsPerTrackVsLUMI->setAxisTitle("scal lumi [10e30 Hz cm^{-2}]", 1);
0544     NumberOfRecHitsPerTrackVsLUMI->setAxisTitle("Mean number of vertices", 2);
0545 
0546     double PVMin = conf->getParameter<double>("PVMin");
0547     double PVMax = conf->getParameter<double>("PVMax");
0548 
0549     histname = "NumberOfGoodPVtxVsLUMI";
0550     NumberOfGoodPVtxVsLUMI = ibooker.bookProfile(histname, histname, LUMIBin, LUMIMin, LUMIMax, PVMin, 3. * PVMax, "");
0551     NumberOfGoodPVtxVsLUMI->setAxisTitle("scal lumi [10e30 Hz cm^{-2}]", 1);
0552     NumberOfGoodPVtxVsLUMI->setAxisTitle("Mean number of vertices", 2);
0553 
0554     histname = "NumberOfGoodPVtxWO0VsLUMI";
0555     NumberOfGoodPVtxWO0VsLUMI =
0556         ibooker.bookProfile(histname, histname, LUMIBin, LUMIMin, LUMIMax, PVMin, 3. * PVMax, "");
0557     NumberOfGoodPVtxWO0VsLUMI->setAxisTitle("scal lumi [10e30 Hz cm^{-2}]", 1);
0558     NumberOfGoodPVtxWO0VsLUMI->setAxisTitle("Mean number of vertices", 2);
0559 
0560     double NClusPxMin = conf->getParameter<double>("NClusPxMin");
0561     double NClusPxMax = conf->getParameter<double>("NClusPxMax");
0562     histname = "NumberOfPixelClustersVsGoodPVtx";
0563     NumberOfPixelClustersVsLUMI =
0564         ibooker.bookProfile(histname, histname, LUMIBin, LUMIMin, LUMIMax, NClusPxMin, 3. * NClusPxMax, "");
0565     NumberOfPixelClustersVsLUMI->setAxisTitle("scal lumi [10e30 Hz cm^{-2}]", 1);
0566     NumberOfPixelClustersVsLUMI->setAxisTitle("Mean number of pixel clusters", 2);
0567 
0568     double NClusStrMin = conf->getParameter<double>("NClusStrMin");
0569     double NClusStrMax = conf->getParameter<double>("NClusStrMax");
0570     histname = "NumberOfStripClustersVsLUMI";
0571     NumberOfStripClustersVsLUMI =
0572         ibooker.bookProfile(histname, histname, LUMIBin, LUMIMin, LUMIMax, NClusStrMin, 3. * NClusStrMax, "");
0573     NumberOfStripClustersVsLUMI->setAxisTitle("scal lumi [10e30 Hz cm^{-2}]", 1);
0574     NumberOfStripClustersVsLUMI->setAxisTitle("Mean number of strip clusters", 2);
0575   }
0576 
0577   if (doPlotsVsBXlumi_) {
0578     ibooker.setCurrentFolder(MEFolderName + "/PUmonitoring");
0579     // get binning from the configuration
0580     edm::ParameterSet BXlumiParameters = conf->getParameter<edm::ParameterSet>("BXlumiSetup");
0581     int BXlumiBin = BXlumiParameters.getParameter<int>("BXlumiBin");
0582     double BXlumiMin = BXlumiParameters.getParameter<double>("BXlumiMin");
0583     double BXlumiMax = BXlumiParameters.getParameter<double>("BXlumiMax");
0584 
0585     histname = "NumberOfTracksVsBXlumi_" + CategoryName;
0586     NumberOfTracksVsBXlumi =
0587         ibooker.bookProfile(histname, histname, BXlumiBin, BXlumiMin, BXlumiMax, TKNoMin, 3. * TKNoMax, "");
0588     NumberOfTracksVsBXlumi->setAxisTitle("lumi BX [10^{30}Hzcm^{-2}]", 1);
0589     NumberOfTracksVsBXlumi->setAxisTitle("Mean number of Tracks", 2);
0590   }
0591 
0592   if (doLumiAnalysis) {
0593     auto scope = DQMStore::IBooker::UseLumiScope(ibooker);
0594     theTrackAnalyzer->initHisto(ibooker, iSetup, *conf);
0595   } else {
0596     theTrackAnalyzer->initHisto(ibooker, iSetup, *conf);
0597   }
0598 
0599   // book the Seed Property histograms
0600   // ---------------------------------------------------------------------------------//
0601 
0602   ibooker.setCurrentFolder(MEFolderName + "/TrackBuilding");
0603 
0604   doSeedNumberPlot = conf->getParameter<bool>("doSeedNumberHisto");
0605   doSeedLumiAnalysis_ = conf->getParameter<bool>("doSeedLumiAnalysis");
0606   doSeedVsClusterPlot = conf->getParameter<bool>("doSeedVsClusterHisto");
0607 
0608   edm::InputTag seedProducer = conf->getParameter<edm::InputTag>("SeedProducer");
0609 
0610   if (doAllSeedPlots || doSeedNumberPlot) {
0611     ibooker.setCurrentFolder(MEFolderName + "/TrackBuilding");
0612     histname = "NumberOfSeeds_" + seedProducer.label() + "_" + CategoryName;
0613     NumberOfSeeds = ibooker.book1D(histname, histname, TKNoSeedBin, TKNoSeedMin, TKNoSeedMax);
0614     NumberOfSeeds->setAxisTitle("Number of Seeds per Event", 1);
0615     NumberOfSeeds->setAxisTitle("Number of Events", 2);
0616 
0617     if (doSeedLumiAnalysis_) {
0618       ibooker.setCurrentFolder(MEFolderName + "/LSanalysis");
0619       auto scope = DQMStore::IBooker::UseLumiScope(ibooker);
0620       histname = "NumberOfSeeds_lumiFlag_" + seedProducer.label() + "_" + CategoryName;
0621       NumberOfSeeds_lumiFlag = ibooker.book1D(histname, histname, TKNoSeedBin, TKNoSeedMin, TKNoSeedMax);
0622       NumberOfSeeds_lumiFlag->setAxisTitle("Number of Seeds per Event", 1);
0623       NumberOfSeeds_lumiFlag->setAxisTitle("Number of Events", 2);
0624     }
0625   }
0626 
0627   if (doAllSeedPlots || doSeedVsClusterPlot) {
0628     ibooker.setCurrentFolder(MEFolderName + "/TrackBuilding");
0629 
0630     ClusterLabels = conf->getParameter<std::vector<std::string> >("ClusterLabels");
0631 
0632     std::vector<double> histoMin, histoMax;
0633     std::vector<int> histoBin;  //these vectors are for max min and nbins in histograms
0634 
0635     int NClusPxBin = conf->getParameter<int>("NClusPxBin");
0636     double NClusPxMin = conf->getParameter<double>("NClusPxMin");
0637     double NClusPxMax = conf->getParameter<double>("NClusPxMax");
0638 
0639     int NClusStrBin = conf->getParameter<int>("NClusStrBin");
0640     double NClusStrMin = conf->getParameter<double>("NClusStrMin");
0641     double NClusStrMax = conf->getParameter<double>("NClusStrMax");
0642 
0643     setMaxMinBin(
0644         histoMin, histoMax, histoBin, NClusStrMin, NClusStrMax, NClusStrBin, NClusPxMin, NClusPxMax, NClusPxBin);
0645 
0646     for (uint i = 0; i < ClusterLabels.size(); i++) {
0647       histname = "SeedsVsClusters_" + seedProducer.label() + "_Vs_" + ClusterLabels[i] + "_" + CategoryName;
0648       SeedsVsClusters.push_back(dynamic_cast<MonitorElement*>(ibooker.book2D(
0649           histname, histname, histoBin[i], histoMin[i], histoMax[i], TKNoSeedBin, TKNoSeedMin, TKNoSeedMax)));
0650       SeedsVsClusters[i]->setAxisTitle("Number of Clusters", 1);
0651       SeedsVsClusters[i]->setAxisTitle("Number of Seeds", 2);
0652       SeedsVsClusters[i]->getTH2F()->SetCanExtend(TH1::kAllAxes);
0653     }
0654   }
0655 
0656   if (doRegionPlots) {
0657     ibooker.setCurrentFolder(MEFolderName + "/TrackBuilding");
0658 
0659     int regionBin = conf->getParameter<int>("RegionSizeBin");
0660     double regionMin = conf->getParameter<double>("RegionSizeMin");
0661     double regionMax = conf->getParameter<double>("RegionSizeMax");
0662 
0663     histname = "TrackingRegionsNumberOf_" + seedProducer.label() + "_" + CategoryName;
0664     NumberOfTrackingRegions = ibooker.book1D(histname, histname, regionBin, regionMin, regionMax);
0665     NumberOfTrackingRegions->setAxisTitle("Number of TrackingRegions per Event", 1);
0666     NumberOfTrackingRegions->setAxisTitle("Number of Events", 2);
0667   }
0668 
0669   if (doTkCandPlots) {
0670     ibooker.setCurrentFolder(MEFolderName + "/TrackBuilding");
0671 
0672     edm::InputTag tcProducer = conf->getParameter<edm::InputTag>("TCProducer");
0673 
0674     histname = "NumberOfTrackCandidates_" + tcProducer.label() + "_" + CategoryName;
0675     NumberOfTrackCandidates = ibooker.book1D(histname, histname, TCNoBin, TCNoMin, TCNoMax);
0676     NumberOfTrackCandidates->setAxisTitle("Number of Track Candidates per Event", 1);
0677     NumberOfTrackCandidates->setAxisTitle("Number of Event", 2);
0678 
0679     histname = "FractionOfCandOverSeeds_" + tcProducer.label() + "_" + CategoryName;
0680     FractionCandidatesOverSeeds = ibooker.book1D(histname, histname, 101, 0., 1.01);
0681     FractionCandidatesOverSeeds->setAxisTitle("Number of Track Candidates / Number of Seeds per Event", 1);
0682     FractionCandidatesOverSeeds->setAxisTitle("Number of Event", 2);
0683   }
0684 
0685   theTrackBuildingAnalyzer->initHisto(ibooker, *conf);
0686 
0687   if (doTrackerSpecific_ || doAllPlots) {
0688     ClusterLabels = conf->getParameter<std::vector<std::string> >("ClusterLabels");
0689 
0690     std::vector<double> histoMin, histoMax;
0691     std::vector<int> histoBin;  //these vectors are for max min and nbins in histograms
0692 
0693     int NClusStrBin = conf->getParameter<int>("NClusStrBin");
0694     double NClusStrMin = conf->getParameter<double>("NClusStrMin");
0695     double NClusStrMax = conf->getParameter<double>("NClusStrMax");
0696 
0697     int NClusPxBin = conf->getParameter<int>("NClusPxBin");
0698     double NClusPxMin = conf->getParameter<double>("NClusPxMin");
0699     double NClusPxMax = conf->getParameter<double>("NClusPxMax");
0700 
0701     edm::ParameterSet ParametersNTrk2D = conf->getParameter<edm::ParameterSet>("NTrk2D");
0702     int NTrk2DBin = ParametersNTrk2D.getParameter<int>("NTrk2DBin");
0703     double NTrk2DMin = ParametersNTrk2D.getParameter<double>("NTrk2DMin");
0704     double NTrk2DMax = ParametersNTrk2D.getParameter<double>("NTrk2DMax");
0705 
0706     setMaxMinBin(
0707         histoMin, histoMax, histoBin, NClusStrMin, NClusStrMax, NClusStrBin, NClusPxMin, NClusPxMax, NClusPxBin);
0708 
0709     ibooker.setCurrentFolder(MEFolderName + "/HitProperties");
0710 
0711     for (uint i = 0; i < ClusterLabels.size(); i++) {
0712       ibooker.setCurrentFolder(MEFolderName + "/HitProperties");
0713       histname = "TracksVs" + ClusterLabels[i] + "Cluster_" + CategoryName;
0714       NumberOfTrkVsClusters.push_back(dynamic_cast<MonitorElement*>(
0715           ibooker.book2D(histname, histname, histoBin[i], histoMin[i], histoMax[i], NTrk2DBin, NTrk2DMin, NTrk2DMax)));
0716       std::string title = "Number of " + ClusterLabels[i] + " Clusters";
0717       if (ClusterLabels[i] == "Tot")
0718         title = "# of Clusters in (Pixel+Strip) Detectors";
0719       NumberOfTrkVsClusters[i]->setAxisTitle(title, 1);
0720       NumberOfTrkVsClusters[i]->setAxisTitle("Number of Tracks", 2);
0721       NumberOfTrkVsClusters[i]->getTH1()->SetCanExtend(TH1::kXaxis);
0722     }
0723   }
0724 
0725   // Initialize the GenericTriggerEventFlag
0726   if (genTriggerEventFlag_->on())
0727     genTriggerEventFlag_->initRun(iRun, iSetup);
0728 }
0729 
0730 /*
0731 // -- BeginRun
0732 //---------------------------------------------------------------------------------//
0733 void TrackingMonitor::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup)
0734 {
0735   
0736 }
0737 */
0738 
0739 // -- Analyse
0740 // ---------------------------------------------------------------------------------//
0741 void TrackingMonitor::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0742   // Filter out events if Trigger Filtering is requested
0743   if (genTriggerEventFlag_->on() && !genTriggerEventFlag_->accept(iEvent, iSetup))
0744     return;
0745   auto const* conf = edm::pset::Registry::instance()->getMapped(confID_);
0746   MEFolderName = conf->getParameter<std::string>("FolderName");
0747   std::string Folder = MEFolderName.substr(0, 2);
0748   float lumi = -1.;
0749   if (forceSCAL_) {
0750     edm::Handle<LumiScalersCollection> lumiScalers = iEvent.getHandle(lumiscalersToken_);
0751     if (lumiScalers.isValid() && !lumiScalers->empty()) {
0752       LumiScalersCollection::const_iterator scalit = lumiScalers->begin();
0753       lumi = scalit->instantLumi();
0754     }
0755   } else {
0756     edm::Handle<OnlineLuminosityRecord> metaData = iEvent.getHandle(metaDataToken_);
0757     if (metaData.isValid())
0758       lumi = metaData->instLumi();
0759   }
0760 
0761   if (doPlotsVsLUMI_ || doAllPlots)
0762     NumberEventsOfVsLUMI->Fill(lumi);
0763 
0764   //  Analyse the tracks
0765   //  if the collection is empty, do not fill anything
0766   // ---------------------------------------------------------------------------------//
0767 
0768   size_t bx = iEvent.bunchCrossing();
0769   if (doPlotsVsBX_ || doAllPlots)
0770     NumberEventsOfVsBX->Fill(bx);
0771 
0772   // get the track collection
0773   edm::Handle<edm::View<reco::Track> > trackHandle = iEvent.getHandle(trackToken_);
0774 
0775   int numberOfTracks_den = 0;
0776   edm::Handle<edm::View<reco::Track> > allTrackHandle = iEvent.getHandle(allTrackToken_);
0777   if (allTrackHandle.isValid()) {
0778     for (edm::View<reco::Track>::const_iterator track = allTrackHandle->begin(); track != allTrackHandle->end();
0779          ++track) {
0780       if (denSelection_(*track))
0781         numberOfTracks_den++;
0782     }
0783   }
0784 
0785   edm::Handle<reco::VertexCollection> pvHandle = iEvent.getHandle(pvSrcToken_);
0786   reco::Vertex const* pv0 = nullptr;
0787   if (pvHandle.isValid()) {
0788     pv0 = &pvHandle->front();
0789     //--- pv fake (the pv collection should have size==1 and the pv==beam spot)
0790     if (pv0->isFake() ||
0791         pv0->tracksSize() == 0
0792         // definition of goodOfflinePrimaryVertex
0793         || pv0->ndof() < pvNDOF_ || pv0->z() > 24.)
0794       pv0 = nullptr;
0795   }
0796 
0797   if (trackHandle.isValid()) {
0798     int numberOfTracks = trackHandle->size();
0799     int numberOfTracks_num = 0;
0800     int numberOfTracks_pv0 = 0;
0801 
0802     const edm::View<reco::Track>& trackCollection = *trackHandle;
0803     // calculate the mean # rechits and layers
0804     int totalRecHits = 0, totalLayers = 0;
0805 
0806     theTrackAnalyzer->setNumberOfGoodVertices(iEvent);
0807     theTrackAnalyzer->setBX(iEvent);
0808     theTrackAnalyzer->setLumi(iEvent, iSetup);
0809     for (edm::View<reco::Track>::const_iterator track = trackCollection.begin(); track != trackCollection.end();
0810          ++track) {
0811       if (doPlotsVsBX_ || doAllPlots)
0812         NumberOfRecHitsPerTrackVsBX->Fill(bx, track->numberOfValidHits());
0813       if (numSelection_(*track)) {
0814         numberOfTracks_num++;
0815         if (pv0 && std::abs(track->dz(pv0->position())) < 0.15)
0816           ++numberOfTracks_pv0;
0817       }
0818 
0819       if (doProfilesVsLS_ || doAllPlots)
0820         NumberOfRecHitsPerTrackVsLS->Fill(static_cast<double>(iEvent.id().luminosityBlock()),
0821                                           track->numberOfValidHits());
0822 
0823       if (doPlotsVsLUMI_ || doAllPlots)
0824         NumberOfRecHitsPerTrackVsLUMI->Fill(lumi, track->numberOfValidHits());
0825 
0826       totalRecHits += track->numberOfValidHits();
0827       totalLayers += track->hitPattern().trackerLayersWithMeasurement();
0828 
0829       // do analysis per track
0830       theTrackAnalyzer->analyze(iEvent, iSetup, *track);
0831     }
0832 
0833     double frac = -1.;
0834     //      if (numberOfAllTracks > 0) frac = static_cast<double>(numberOfTracks)/static_cast<double>(numberOfAllTracks);
0835     if (numberOfTracks_den > 0)
0836       frac = static_cast<double>(numberOfTracks_num) / static_cast<double>(numberOfTracks_den);
0837 
0838     if (doGeneralPropertiesPlots_ || doAllPlots) {
0839       NumberOfTracks->Fill(double(numberOfTracks));
0840 
0841       if (Folder == "Tr") {
0842         NumberofTracks_Hardvtx->Fill(double(numberOfTracks_pv0));
0843         NumberOfTracks_PUvtx->Fill(double(numberOfTracks - numberOfTracks_pv0));
0844         NumberofTracks_Hardvtx_PUvtx->Fill(0.5, double(numberOfTracks_pv0));
0845         NumberofTracks_Hardvtx_PUvtx->Fill(1.5, double(numberOfTracks - numberOfTracks_pv0));
0846       }
0847 
0848       if (doPlotsVsBX_ || doAllPlots)
0849         NumberOfTracksVsBX->Fill(bx, numberOfTracks);
0850       if (doPlotsVsLUMI_ || doAllPlots)
0851         NumberOfTracksVsLUMI->Fill(lumi, numberOfTracks);
0852       if (doFractionPlot_) {
0853         FractionOfGoodTracks->Fill(frac);
0854 
0855         if (doFractionPlot_) {
0856           if (doPlotsVsBX_ || doAllPlots)
0857             GoodTracksFractionVsBX->Fill(bx, frac);
0858           if (doPlotsVsLUMI_ || doAllPlots)
0859             GoodTracksFractionVsLUMI->Fill(lumi, frac);
0860         }
0861       }
0862       if (numberOfTracks > 0) {
0863         double meanRecHits = static_cast<double>(totalRecHits) / static_cast<double>(numberOfTracks);
0864         double meanLayers = static_cast<double>(totalLayers) / static_cast<double>(numberOfTracks);
0865         NumberOfMeanRecHitsPerTrack->Fill(meanRecHits);
0866         NumberOfMeanLayersPerTrack->Fill(meanLayers);
0867       }
0868     }
0869 
0870     if (doProfilesVsLS_ || doAllPlots) {
0871       float nLS = static_cast<double>(iEvent.id().luminosityBlock());
0872       NumberEventsOfVsLS->Fill(nLS);
0873       NumberOfTracksVsLS->Fill(nLS, numberOfTracks);
0874       if (doFractionPlot_)
0875         GoodTracksFractionVsLS->Fill(nLS, frac);
0876     }
0877 
0878     if (doLumiAnalysis) {
0879       NumberOfTracks_lumiFlag->Fill(numberOfTracks);
0880     }
0881 
0882     //  Analyse the Track Building variables
0883     //  if the collection is empty, do not fill anything
0884     // ---------------------------------------------------------------------------------//
0885 
0886     // fill the TrackCandidate info
0887     if (doTkCandPlots) {
0888       // magnetic field
0889       MagneticField const& theMF = iSetup.getData(magneticFieldToken_);
0890 
0891       // get the candidate collection
0892       edm::Handle<TrackCandidateCollection> theTCHandle = iEvent.getHandle(trackCandidateToken_);
0893       const TrackCandidateCollection& theTCCollection = *theTCHandle;
0894 
0895       if (theTCHandle.isValid()) {
0896         // get the beam spot
0897         edm::Handle<reco::BeamSpot> recoBeamSpotHandle = iEvent.getHandle(bsSrcToken_);
0898         const reco::BeamSpot& bs = *recoBeamSpotHandle;
0899 
0900         NumberOfTrackCandidates->Fill(theTCCollection.size());
0901 
0902         // get the seed collection
0903         edm::Handle<edm::View<TrajectorySeed> > seedHandle = iEvent.getHandle(seedToken_);
0904         const edm::View<TrajectorySeed>& seedCollection = *seedHandle;
0905         if (seedHandle.isValid() && !seedCollection.empty())
0906           FractionCandidatesOverSeeds->Fill(double(theTCCollection.size()) / double(seedCollection.size()));
0907 
0908         TransientTrackingRecHitBuilder const& theTTRHBuilder = iSetup.getData(transientTrackingRecHitBuilderToken_);
0909         for (TrackCandidateCollection::const_iterator cand = theTCCollection.begin(); cand != theTCCollection.end();
0910              ++cand) {
0911           theTrackBuildingAnalyzer->analyze(iEvent, iSetup, *cand, bs, theMF, theTTRHBuilder);
0912         }
0913       } else {
0914         edm::LogWarning("TrackingMonitor") << "No Track Candidates in the event.  Not filling associated histograms";
0915       }
0916 
0917       if (doMVAPlots) {
0918         // Get MVA and quality mask collections
0919         std::vector<const MVACollection*> mvaCollections;
0920         std::vector<const QualityMaskCollection*> qualityMaskCollections;
0921 
0922         edm::Handle<edm::View<reco::Track> > htracks = iEvent.getHandle(mvaTrackToken_);
0923 
0924         edm::Handle<MVACollection> hmva;
0925         edm::Handle<QualityMaskCollection> hqual;
0926         for (const auto& tokenTpl : mvaQualityTokens_) {
0927           hmva = iEvent.getHandle(std::get<0>(tokenTpl));
0928           hqual = iEvent.getHandle(std::get<1>(tokenTpl));
0929           mvaCollections.push_back(hmva.product());
0930           qualityMaskCollections.push_back(hqual.product());
0931         }
0932         theTrackBuildingAnalyzer->analyze(*htracks, mvaCollections, qualityMaskCollections);
0933       }
0934     }
0935 
0936     //plots for trajectory seeds
0937 
0938     if (doAllSeedPlots || doSeedNumberPlot || doSeedVsClusterPlot || runTrackBuildingAnalyzerForSeed) {
0939       // get the seed collection
0940       edm::Handle<edm::View<TrajectorySeed> > seedHandle = iEvent.getHandle(seedToken_);
0941 
0942       // fill the seed info
0943       if (seedHandle.isValid()) {
0944         const auto& seedCollection = *seedHandle;
0945 
0946         if (doAllSeedPlots || doSeedNumberPlot) {
0947           NumberOfSeeds->Fill(seedCollection.size());
0948           if (doSeedLumiAnalysis_)
0949             NumberOfSeeds_lumiFlag->Fill(seedCollection.size());
0950         }
0951 
0952         if (doAllSeedPlots || doSeedVsClusterPlot) {
0953           std::vector<int> NClus;
0954           setNclus(iEvent, NClus);
0955           for (uint i = 0; i < ClusterLabels.size(); i++) {
0956             SeedsVsClusters[i]->Fill(NClus[i], seedCollection.size());
0957           }
0958         }
0959 
0960         if (doAllSeedPlots || runTrackBuildingAnalyzerForSeed) {
0961           edm::Handle<std::vector<SeedStopInfo> > stopHandle = iEvent.getHandle(seedStopInfoToken_);
0962           const auto& seedStopInfo = *stopHandle;
0963 
0964           if (seedStopInfo.size() == seedCollection.size()) {
0965             //here duplication of mag field and be informations is needed to allow seed and track cand histos to be independent
0966             // magnetic field
0967             MagneticField const& theMF = iSetup.getData(magneticFieldToken_);
0968 
0969             // get the beam spot
0970             edm::Handle<reco::BeamSpot> recoBeamSpotHandle = iEvent.getHandle(bsSrcToken_);
0971             const reco::BeamSpot& bs = *recoBeamSpotHandle;
0972 
0973             TransientTrackingRecHitBuilder const& theTTRHBuilder = iSetup.getData(transientTrackingRecHitBuilderToken_);
0974             for (size_t i = 0; i < seedCollection.size(); ++i) {
0975               theTrackBuildingAnalyzer->analyze(
0976                   iEvent, iSetup, seedCollection[i], seedStopInfo[i], bs, theMF, theTTRHBuilder);
0977             }
0978           } else {
0979             edm::LogWarning("TrackingMonitor")
0980                 << "Seed collection size (" << seedCollection.size()
0981                 << ") differs from seed stop info collection size (" << seedStopInfo.size()
0982                 << "). This is a sign of inconsistency in the configuration. Not filling associated histograms.";
0983           }
0984         }
0985 
0986       } else {
0987         edm::LogWarning("TrackingMonitor") << "No Trajectory seeds in the event.  Not filling associated histograms";
0988       }
0989     }
0990 
0991     // plots for tracking regions
0992     if (doRegionPlots) {
0993       if (!regionToken_.isUninitialized()) {
0994         edm::Handle<edm::OwnVector<TrackingRegion> > hregions = iEvent.getHandle(regionToken_);
0995         const auto& regions = *hregions;
0996         NumberOfTrackingRegions->Fill(regions.size());
0997 
0998         theTrackBuildingAnalyzer->analyze(regions);
0999       } else if (!regionLayerSetsToken_.isUninitialized()) {
1000         edm::Handle<TrackingRegionsSeedingLayerSets> hregions = iEvent.getHandle(regionLayerSetsToken_);
1001         const auto& regions = *hregions;
1002         NumberOfTrackingRegions->Fill(regions.regionsSize());
1003 
1004         theTrackBuildingAnalyzer->analyze(regions);
1005       }
1006 
1007       if (doRegionCandidatePlots) {
1008         edm::Handle<reco::CandidateView> hcandidates = iEvent.getHandle(regionCandidateToken_);
1009         theTrackBuildingAnalyzer->analyze(*hcandidates);
1010       }
1011     }
1012 
1013     if (doTrackerSpecific_ || doAllPlots) {
1014       std::vector<int> NClus;
1015       setNclus(iEvent, NClus);
1016       for (uint i = 0; i < ClusterLabels.size(); i++) {
1017         NumberOfTrkVsClusters[i]->Fill(NClus[i], numberOfTracks);
1018       }
1019     }
1020 
1021     if (doPUmonitoring_) {
1022       // do vertex monitoring
1023       for (size_t i = 0; i < theVertexMonitor.size(); i++)
1024         theVertexMonitor[i]->analyze(iEvent, iSetup);
1025     }
1026     if (doPlotsVsGoodPVtx_) {
1027       size_t totalNumGoodPV = 0;
1028       if (pvHandle.isValid()) {
1029         for (reco::VertexCollection::const_iterator pv = pvHandle->begin(); pv != pvHandle->end(); ++pv) {
1030           //--- pv fake (the pv collection should have size==1 and the pv==beam spot)
1031           if (pv->isFake() || pv->tracksSize() == 0)
1032             continue;
1033 
1034           // definition of goodOfflinePrimaryVertex
1035           if (pv->ndof() < pvNDOF_ || pv->z() > 24.)
1036             continue;
1037           totalNumGoodPV++;
1038         }
1039 
1040         NumberEventsOfVsGoodPVtx->Fill(float(totalNumGoodPV));
1041         NumberOfTracksVsGoodPVtx->Fill(float(totalNumGoodPV), numberOfTracks);
1042         if (totalNumGoodPV > 1)
1043           NumberOfTracksVsPUPVtx->Fill(totalNumGoodPV - 1,
1044                                        double(numberOfTracks - numberOfTracks_pv0) / double(totalNumGoodPV - 1));
1045         NumberOfPVtxVsGoodPVtx->Fill(float(totalNumGoodPV), pvHandle->size());
1046 
1047         for (edm::View<reco::Track>::const_iterator track = trackCollection.begin(); track != trackCollection.end();
1048              ++track) {
1049           NumberOfRecHitsPerTrackVsGoodPVtx->Fill(float(totalNumGoodPV), track->numberOfValidHits());
1050         }
1051 
1052         if (doProfilesVsLS_ || doAllPlots)
1053           NumberOfGoodPVtxVsLS->Fill(static_cast<double>(iEvent.id().luminosityBlock()), totalNumGoodPV);
1054         if (doPlotsVsBX_ || doAllPlots)
1055           NumberOfGoodPVtxVsBX->Fill(bx, float(totalNumGoodPV));
1056 
1057         if (doFractionPlot_)
1058           GoodTracksFractionVsGoodPVtx->Fill(float(totalNumGoodPV), frac);
1059 
1060         if (doPlotsVsLUMI_ || doAllPlots)
1061           NumberOfGoodPVtxVsLUMI->Fill(lumi, float(totalNumGoodPV));
1062       }
1063 
1064       std::vector<int> NClus;
1065       setNclus(iEvent, NClus);
1066       std::ostringstream ss;
1067       ss << "VI stat " << totalNumGoodPV << ' ' << numberOfTracks;
1068       for (uint i = 0; i < ClusterLabels.size(); i++) {
1069         ss << ' ' << NClus[i];
1070         if (doPlotsVsLUMI_ || doAllPlots) {
1071           if (ClusterLabels[i] == "Pix")
1072             NumberOfPixelClustersVsLUMI->Fill(lumi, NClus[i]);
1073           if (ClusterLabels[i] == "Strip")
1074             NumberOfStripClustersVsLUMI->Fill(lumi, NClus[i]);
1075         }
1076         if (ClusterLabels[i] == "Pix")
1077           NumberOfPixelClustersVsGoodPVtx->Fill(float(totalNumGoodPV), NClus[i]);
1078         if (ClusterLabels[i] == "Strip")
1079           NumberOfStripClustersVsGoodPVtx->Fill(float(totalNumGoodPV), NClus[i]);
1080       }
1081       COUT(MEFolderName) << ss.str() << std::endl;
1082       if (doPlotsVsBXlumi_) {
1083         double bxlumi = theLumiDetails_->getValue(iEvent);
1084         NumberOfTracksVsBXlumi->Fill(bxlumi, numberOfTracks);
1085       }
1086 
1087       if (doProfilesVsLS_ || doAllPlots)
1088         if (totalNumGoodPV != 0)
1089           NumberOfGoodPVtxWO0VsLS->Fill(static_cast<double>(iEvent.id().luminosityBlock()), float(totalNumGoodPV));
1090       if (doPlotsVsBX_ || doAllPlots)
1091         if (totalNumGoodPV != 0)
1092           NumberOfGoodPVtxWO0VsBX->Fill(bx, float(totalNumGoodPV));
1093       if (doPlotsVsLUMI_ || doAllPlots)
1094         if (totalNumGoodPV != 0)
1095           NumberOfGoodPVtxWO0VsLUMI->Fill(lumi, float(totalNumGoodPV));
1096 
1097     }  // PU monitoring
1098 
1099   }  // trackHandle is valid
1100 }
1101 
1102 void TrackingMonitor::setMaxMinBin(std::vector<double>& arrayMin,
1103                                    std::vector<double>& arrayMax,
1104                                    std::vector<int>& arrayBin,
1105                                    double smin,
1106                                    double smax,
1107                                    int sbin,
1108                                    double pmin,
1109                                    double pmax,
1110                                    int pbin) {
1111   arrayMin.resize(ClusterLabels.size());
1112   arrayMax.resize(ClusterLabels.size());
1113   arrayBin.resize(ClusterLabels.size());
1114 
1115   for (uint i = 0; i < ClusterLabels.size(); ++i) {
1116     if (ClusterLabels[i] == "Pix") {
1117       arrayMin[i] = pmin;
1118       arrayMax[i] = pmax;
1119       arrayBin[i] = pbin;
1120     } else if (ClusterLabels[i] == "Strip") {
1121       arrayMin[i] = smin;
1122       arrayMax[i] = smax;
1123       arrayBin[i] = sbin;
1124     } else if (ClusterLabels[i] == "Tot") {
1125       arrayMin[i] = smin;
1126       arrayMax[i] = smax + pmax;
1127       arrayBin[i] = sbin;
1128     } else {
1129       edm::LogWarning("TrackingMonitor") << "Cluster Label " << ClusterLabels[i]
1130                                          << " not defined, using strip parameters ";
1131       arrayMin[i] = smin;
1132       arrayMax[i] = smax;
1133       arrayBin[i] = sbin;
1134     }
1135   }
1136 }
1137 
1138 void TrackingMonitor::setNclus(const edm::Event& iEvent, std::vector<int>& arrayNclus) {
1139   int ncluster_pix = -1;
1140   int ncluster_strip = -1;
1141 
1142   edm::Handle<edmNew::DetSetVector<SiStripCluster> > strip_clusters = iEvent.getHandle(stripClustersToken_);
1143   edm::Handle<edmNew::DetSetVector<SiPixelCluster> > pixel_clusters = iEvent.getHandle(pixelClustersToken_);
1144 
1145   if (strip_clusters.isValid() && pixel_clusters.isValid()) {
1146     ncluster_pix = (*pixel_clusters).dataSize();
1147     ncluster_strip = (*strip_clusters).dataSize();
1148   }
1149 
1150   arrayNclus.resize(ClusterLabels.size());
1151   for (uint i = 0; i < ClusterLabels.size(); ++i) {
1152     if (ClusterLabels[i] == "Pix")
1153       arrayNclus[i] = ncluster_pix;
1154     else if (ClusterLabels[i] == "Strip")
1155       arrayNclus[i] = ncluster_strip;
1156     else if (ClusterLabels[i] == "Tot")
1157       arrayNclus[i] = ncluster_pix + ncluster_strip;
1158     else {
1159       edm::LogWarning("TrackingMonitor") << "Cluster Label " << ClusterLabels[i]
1160                                          << " not defined using stri parametrs ";
1161       arrayNclus[i] = ncluster_strip;
1162     }
1163   }
1164 }
1165 DEFINE_FWK_MODULE(TrackingMonitor);