Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2023-06-02 00:51:03 

0001 #include "Validation/HGCalValidation/interface/HGCalValidator.h"
0002 
0003 #include "SimCalorimetry/HGCalAssociatorProducers/interface/AssociatorTools.h"
0004 
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 
0008 using namespace std;
0009 using namespace edm;
0010 
0011 HGCalValidator::HGCalValidator(const edm::ParameterSet& pset)
0012     : caloGeomToken_(esConsumes<CaloGeometry, CaloGeometryRecord>()),
0013       label_lcl(pset.getParameter<edm::InputTag>("label_lcl")),
0014       label_tst(pset.getParameter<std::vector<edm::InputTag>>("label_tst")),
0015       label_simTS(pset.getParameter<edm::InputTag>("label_simTS")),
0016       label_simTSFromCP(pset.getParameter<edm::InputTag>("label_simTSFromCP")),
0017       associator_(pset.getUntrackedParameter<edm::InputTag>("associator")),
0018       associatorSim_(pset.getUntrackedParameter<edm::InputTag>("associatorSim")),
0019       SaveGeneralInfo_(pset.getUntrackedParameter<bool>("SaveGeneralInfo")),
0020       doCaloParticlePlots_(pset.getUntrackedParameter<bool>("doCaloParticlePlots")),
0021       doCaloParticleSelection_(pset.getUntrackedParameter<bool>("doCaloParticleSelection")),
0022       doSimClustersPlots_(pset.getUntrackedParameter<bool>("doSimClustersPlots")),
0023       label_SimClustersPlots_(pset.getParameter<edm::InputTag>("label_SimClusters")),
0024       label_SimClustersLevel_(pset.getParameter<edm::InputTag>("label_SimClustersLevel")),
0025       doLayerClustersPlots_(pset.getUntrackedParameter<bool>("doLayerClustersPlots")),
0026       label_layerClustersPlots_(pset.getParameter<edm::InputTag>("label_layerClusterPlots")),
0027       label_LCToCPLinking_(pset.getParameter<edm::InputTag>("label_LCToCPLinking")),
0028       doTrackstersPlots_(pset.getUntrackedParameter<bool>("doTrackstersPlots")),
0029       label_TS_(pset.getParameter<std::string>("label_TS")),
0030       label_TSToCPLinking_(pset.getParameter<std::string>("label_TSToCPLinking")),
0031       label_TSToSTSPR_(pset.getParameter<std::string>("label_TSToSTSPR")),
0032       label_clustersmask(pset.getParameter<std::vector<edm::InputTag>>("LayerClustersInputMask")),
0033       cummatbudinxo_(pset.getParameter<edm::FileInPath>("cummatbudinxo")) {
0034   //In this way we can easily generalize to associations between other objects also.
0035   const edm::InputTag& label_cp_effic_tag = pset.getParameter<edm::InputTag>("label_cp_effic");
0036   const edm::InputTag& label_cp_fake_tag = pset.getParameter<edm::InputTag>("label_cp_fake");
0037 
0038   label_cp_effic = consumes<std::vector<CaloParticle>>(label_cp_effic_tag);
0039   label_cp_fake = consumes<std::vector<CaloParticle>>(label_cp_fake_tag);
0040 
0041   simVertices_ = consumes<std::vector<SimVertex>>(pset.getParameter<edm::InputTag>("simVertices"));
0042 
0043   for (auto& itag : label_clustersmask) {
0044     clustersMaskTokens_.push_back(consumes<std::vector<float>>(itag));
0045   }
0046 
0047   associatorMapSimtR = consumes<hgcal::SimToRecoCollectionWithSimClusters>(associatorSim_);
0048   associatorMapRtSim = consumes<hgcal::RecoToSimCollectionWithSimClusters>(associatorSim_);
0049 
0050   simTrackstersMap_ = consumes<std::map<uint, std::vector<uint>>>(edm::InputTag("ticlSimTracksters"));
0051 
0052   hitMap_ = consumes<std::unordered_map<DetId, const HGCRecHit*>>(edm::InputTag("hgcalRecHitMapProducer"));
0053 
0054   simClusters_ = consumes<std::vector<SimCluster>>(pset.getParameter<edm::InputTag>("label_scl"));
0055 
0056   layerclusters_ = consumes<reco::CaloClusterCollection>(label_lcl);
0057 
0058   for (auto& itag : label_tst) {
0059     label_tstTokens.push_back(consumes<ticl::TracksterCollection>(itag));
0060   }
0061 
0062   simTracksters_ = consumes<ticl::TracksterCollection>(label_simTS);
0063   simTracksters_fromCPs_ = consumes<ticl::TracksterCollection>(label_simTSFromCP);
0064 
0065   associatorMapRtS = consumes<hgcal::RecoToSimCollection>(associator_);
0066   associatorMapStR = consumes<hgcal::SimToRecoCollection>(associator_);
0067 
0068   cpSelector = CaloParticleSelector(pset.getParameter<double>("ptMinCP"),
0069                                     pset.getParameter<double>("ptMaxCP"),
0070                                     pset.getParameter<double>("minRapidityCP"),
0071                                     pset.getParameter<double>("maxRapidityCP"),
0072                                     pset.getParameter<double>("lipCP"),
0073                                     pset.getParameter<double>("tipCP"),
0074                                     pset.getParameter<int>("minHitCP"),
0075                                     pset.getParameter<int>("maxSimClustersCP"),
0076                                     pset.getParameter<bool>("signalOnlyCP"),
0077                                     pset.getParameter<bool>("intimeOnlyCP"),
0078                                     pset.getParameter<bool>("chargedOnlyCP"),
0079                                     pset.getParameter<bool>("stableOnlyCP"),
0080                                     pset.getParameter<bool>("notConvertedOnlyCP"),
0081                                     pset.getParameter<std::vector<int>>("pdgIdCP"));
0082 
0083   tools_.reset(new hgcal::RecHitTools());
0084 
0085   particles_to_monitor_ = pset.getParameter<std::vector<int>>("pdgIdCP");
0086   totallayers_to_monitor_ = pset.getParameter<int>("totallayers_to_monitor");
0087   thicknesses_to_monitor_ = pset.getParameter<std::vector<int>>("thicknesses_to_monitor");
0088 
0089   //For the material budget file here
0090   std::ifstream fmb(cummatbudinxo_.fullPath().c_str());
0091   double thelay = 0.;
0092   double mbg = 0.;
0093   for (unsigned ilayer = 1; ilayer <= totallayers_to_monitor_; ++ilayer) {
0094     fmb >> thelay >> mbg;
0095     cummatbudg.insert(std::pair<double, double>(thelay, mbg));
0096   }
0097 
0098   fmb.close();
0099 
0100   ParameterSet psetForHistoProducerAlgo = pset.getParameter<ParameterSet>("histoProducerAlgoBlock");
0101   histoProducerAlgo_ = std::make_unique<HGVHistoProducerAlgo>(psetForHistoProducerAlgo);
0102 
0103   dirName_ = pset.getParameter<std::string>("dirName");
0104 }
0105 
0106 HGCalValidator::~HGCalValidator() {}
0107 
0108 void HGCalValidator::bookHistograms(DQMStore::IBooker& ibook,
0109                                     edm::Run const&,
0110                                     edm::EventSetup const& setup,
0111                                     Histograms& histograms) const {
0112   if (SaveGeneralInfo_) {
0113     ibook.cd();
0114     ibook.setCurrentFolder(dirName_ + "GeneralInfo");
0115     histoProducerAlgo_->bookInfo(ibook, histograms.histoProducerAlgo);
0116   }
0117 
0118   if (doCaloParticlePlots_) {
0119     ibook.cd();
0120 
0121     for (auto const particle : particles_to_monitor_) {
0122       ibook.setCurrentFolder(dirName_ + "SelectedCaloParticles/" + std::to_string(particle));
0123       histoProducerAlgo_->bookCaloParticleHistos(
0124           ibook, histograms.histoProducerAlgo, particle, totallayers_to_monitor_);
0125     }
0126     ibook.cd();
0127     ibook.setCurrentFolder(dirName_);
0128   }
0129 
0130   //Booking histograms concerning with simClusters
0131   if (doSimClustersPlots_) {
0132     ibook.cd();
0133     ibook.setCurrentFolder(dirName_ + label_SimClustersPlots_.label() + "/" + label_SimClustersLevel_.label());
0134     histoProducerAlgo_->bookSimClusterHistos(
0135         ibook, histograms.histoProducerAlgo, totallayers_to_monitor_, thicknesses_to_monitor_);
0136 
0137     for (unsigned int ws = 0; ws < label_clustersmask.size(); ws++) {
0138       ibook.cd();
0139       InputTag algo = label_clustersmask[ws];
0140       string dirName = dirName_ + label_SimClustersPlots_.label() + "/";
0141       if (!algo.process().empty())
0142         dirName += algo.process() + "_";
0143       LogDebug("HGCalValidator") << dirName << "\n";
0144       if (!algo.label().empty())
0145         dirName += algo.label() + "_";
0146       LogDebug("HGCalValidator") << dirName << "\n";
0147       if (!algo.instance().empty())
0148         dirName += algo.instance() + "_";
0149       LogDebug("HGCalValidator") << dirName << "\n";
0150 
0151       if (!dirName.empty()) {
0152         dirName.resize(dirName.size() - 1);
0153       }
0154 
0155       LogDebug("HGCalValidator") << dirName << "\n";
0156 
0157       ibook.setCurrentFolder(dirName);
0158 
0159       histoProducerAlgo_->bookSimClusterAssociationHistos(
0160           ibook, histograms.histoProducerAlgo, totallayers_to_monitor_, thicknesses_to_monitor_);
0161     }  //end of loop over masks
0162   }    //if for simCluster plots
0163 
0164   //Booking histograms concerning with hgcal layer clusters
0165   if (doLayerClustersPlots_) {
0166     ibook.cd();
0167     ibook.setCurrentFolder(dirName_ + label_layerClustersPlots_.label() + "/ClusterLevel");
0168     histoProducerAlgo_->bookClusterHistos_ClusterLevel(ibook,
0169                                                        histograms.histoProducerAlgo,
0170                                                        totallayers_to_monitor_,
0171                                                        thicknesses_to_monitor_,
0172                                                        cummatbudinxo_.fullPath());
0173     ibook.cd();
0174     ibook.setCurrentFolder(dirName_ + label_layerClustersPlots_.label() + "/" + label_LCToCPLinking_.label());
0175     histoProducerAlgo_->bookClusterHistos_LCtoCP_association(
0176         ibook, histograms.histoProducerAlgo, totallayers_to_monitor_);
0177 
0178     ibook.cd();
0179     ibook.setCurrentFolder(dirName_ + label_layerClustersPlots_.label() + "/CellLevel");
0180     histoProducerAlgo_->bookClusterHistos_CellLevel(
0181         ibook, histograms.histoProducerAlgo, totallayers_to_monitor_, thicknesses_to_monitor_);
0182   }
0183 
0184   //Booking histograms for Tracksters
0185   for (unsigned int www = 0; www < label_tst.size(); www++) {
0186     ibook.cd();
0187     InputTag algo = label_tst[www];
0188     string dirName = dirName_;
0189     if (!algo.process().empty())
0190       dirName += algo.process() + "_";
0191     LogDebug("HGCalValidator") << dirName << "\n";
0192     if (!algo.label().empty())
0193       dirName += algo.label() + "_";
0194     LogDebug("HGCalValidator") << dirName << "\n";
0195     if (!algo.instance().empty())
0196       dirName += algo.instance() + "_";
0197     LogDebug("HGCalValidator") << dirName << "\n";
0198 
0199     if (!dirName.empty()) {
0200       dirName.resize(dirName.size() - 1);
0201     }
0202 
0203     LogDebug("HGCalValidator") << dirName << "\n";
0204 
0205     ibook.setCurrentFolder(dirName);
0206 
0207     // Booking histograms concerning HGCal tracksters
0208     if (doTrackstersPlots_) {
0209       // Generic histos
0210       ibook.setCurrentFolder(dirName + "/" + label_TS_);
0211       histoProducerAlgo_->bookTracksterHistos(ibook, histograms.histoProducerAlgo, totallayers_to_monitor_);
0212       // CP Linking
0213       ibook.setCurrentFolder(dirName + "/" + label_TSToCPLinking_);
0214       histoProducerAlgo_->bookTracksterSTSHistos(
0215           ibook, histograms.histoProducerAlgo, HGVHistoProducerAlgo::validationType::Linking);
0216       // SimTracksters Pattern Recognition
0217       ibook.setCurrentFolder(dirName + "/" + label_TSToSTSPR_);
0218       histoProducerAlgo_->bookTracksterSTSHistos(
0219           ibook, histograms.histoProducerAlgo, HGVHistoProducerAlgo::validationType::PatternRecognition);
0220     }
0221   }  //end of booking Tracksters loop
0222 }
0223 
0224 void HGCalValidator::cpParametersAndSelection(const Histograms& histograms,
0225                                               std::vector<CaloParticle> const& cPeff,
0226                                               std::vector<SimVertex> const& simVertices,
0227                                               std::vector<size_t>& selected_cPeff,
0228                                               unsigned int layers,
0229                                               std::unordered_map<DetId, const HGCRecHit*> const& hitMap) const {
0230   selected_cPeff.reserve(cPeff.size());
0231 
0232   size_t j = 0;
0233   for (auto const& caloParticle : cPeff) {
0234     int id = caloParticle.pdgId();
0235 
0236     if (!doCaloParticleSelection_ || (doCaloParticleSelection_ && cpSelector(caloParticle, simVertices))) {
0237       selected_cPeff.push_back(j);
0238       if (doCaloParticlePlots_) {
0239         histoProducerAlgo_->fill_caloparticle_histos(
0240             histograms.histoProducerAlgo, id, caloParticle, simVertices, layers, hitMap);
0241       }
0242     }
0243     ++j;
0244   }  //end of loop over caloparticles
0245 }
0246 
0247 void HGCalValidator::dqmAnalyze(const edm::Event& event,
0248                                 const edm::EventSetup& setup,
0249                                 const Histograms& histograms) const {
0250   using namespace reco;
0251 
0252   LogDebug("HGCalValidator") << "\n===================================================="
0253                              << "\n"
0254                              << "Analyzing new event"
0255                              << "\n"
0256                              << "====================================================\n"
0257                              << "\n";
0258 
0259   edm::Handle<std::vector<SimVertex>> simVerticesHandle;
0260   event.getByToken(simVertices_, simVerticesHandle);
0261   std::vector<SimVertex> const& simVertices = *simVerticesHandle;
0262 
0263   edm::Handle<std::vector<CaloParticle>> caloParticleHandle;
0264   event.getByToken(label_cp_effic, caloParticleHandle);
0265   std::vector<CaloParticle> const& caloParticles = *caloParticleHandle;
0266 
0267   edm::Handle<ticl::TracksterCollection> simTracksterHandle;
0268   event.getByToken(simTracksters_, simTracksterHandle);
0269   ticl::TracksterCollection const& simTracksters = *simTracksterHandle;
0270 
0271   edm::Handle<ticl::TracksterCollection> simTracksterFromCPHandle;
0272   event.getByToken(simTracksters_fromCPs_, simTracksterFromCPHandle);
0273   ticl::TracksterCollection const& simTrackstersFromCPs = *simTracksterFromCPHandle;
0274 
0275   edm::Handle<std::map<uint, std::vector<uint>>> simTrackstersMapHandle;
0276   event.getByToken(simTrackstersMap_, simTrackstersMapHandle);
0277   const std::map<uint, std::vector<uint>> cpToSc_SimTrackstersMap = *simTrackstersMapHandle;
0278 
0279   edm::ESHandle<CaloGeometry> geom = setup.getHandle(caloGeomToken_);
0280   tools_->setGeometry(*geom);
0281   histoProducerAlgo_->setRecHitTools(tools_);
0282 
0283   edm::Handle<hgcal::SimToRecoCollection> simtorecoCollectionH;
0284   event.getByToken(associatorMapStR, simtorecoCollectionH);
0285   auto simRecColl = *simtorecoCollectionH;
0286   edm::Handle<hgcal::RecoToSimCollection> recotosimCollectionH;
0287   event.getByToken(associatorMapRtS, recotosimCollectionH);
0288   auto recSimColl = *recotosimCollectionH;
0289 
0290   edm::Handle<std::unordered_map<DetId, const HGCRecHit*>> hitMapHandle;
0291   event.getByToken(hitMap_, hitMapHandle);
0292   const std::unordered_map<DetId, const HGCRecHit*>* hitMap = &*hitMapHandle;
0293 
0294   //Some general info on layers etc.
0295   if (SaveGeneralInfo_) {
0296     histoProducerAlgo_->fill_info_histos(histograms.histoProducerAlgo, totallayers_to_monitor_);
0297   }
0298 
0299   std::vector<size_t> cPIndices;
0300   //Consider CaloParticles coming from the hard scatterer
0301   //excluding the PU contribution and save the indices.
0302   removeCPFromPU(caloParticles, cPIndices);
0303 
0304   // ##############################################
0305   // Fill caloparticles histograms
0306   // ##############################################
0307   // HGCRecHit are given to select the SimHits which are also reconstructed
0308   LogTrace("HGCalValidator") << "\n# of CaloParticles: " << caloParticles.size() << "\n" << std::endl;
0309   std::vector<size_t> selected_cPeff;
0310   cpParametersAndSelection(histograms, caloParticles, simVertices, selected_cPeff, totallayers_to_monitor_, *hitMap);
0311 
0312   //get collections from the event
0313   //simClusters
0314   edm::Handle<std::vector<SimCluster>> simClustersHandle;
0315   event.getByToken(simClusters_, simClustersHandle);
0316   std::vector<SimCluster> const& simClusters = *simClustersHandle;
0317 
0318   //Layer clusters
0319   edm::Handle<reco::CaloClusterCollection> clusterHandle;
0320   event.getByToken(layerclusters_, clusterHandle);
0321   const reco::CaloClusterCollection& clusters = *clusterHandle;
0322 
0323   auto nSimClusters = simClusters.size();
0324   std::vector<size_t> sCIndices;
0325   //There shouldn't be any SimTracks from different crossings, but maybe they will be added later.
0326   //At the moment there should be one SimTrack in each SimCluster.
0327   for (unsigned int scId = 0; scId < nSimClusters; ++scId) {
0328     if (simClusters[scId].g4Tracks()[0].eventId().event() != 0 or
0329         simClusters[scId].g4Tracks()[0].eventId().bunchCrossing() != 0) {
0330       LogDebug("HGCalValidator") << "Excluding SimClusters from event: "
0331                                  << simClusters[scId].g4Tracks()[0].eventId().event()
0332                                  << " with BX: " << simClusters[scId].g4Tracks()[0].eventId().bunchCrossing()
0333                                  << std::endl;
0334       continue;
0335     }
0336     sCIndices.emplace_back(scId);
0337   }
0338 
0339   // ##############################################
0340   // Fill simCluster histograms
0341   // ##############################################
0342   if (doSimClustersPlots_) {
0343     histoProducerAlgo_->fill_simCluster_histos(
0344         histograms.histoProducerAlgo, simClusters, totallayers_to_monitor_, thicknesses_to_monitor_);
0345 
0346     for (unsigned int ws = 0; ws < label_clustersmask.size(); ws++) {
0347       const auto& inputClusterMask = event.get(clustersMaskTokens_[ws]);
0348 
0349       edm::Handle<hgcal::SimToRecoCollectionWithSimClusters> simtorecoCollectionH;
0350       event.getByToken(associatorMapSimtR, simtorecoCollectionH);
0351       auto simRecColl = *simtorecoCollectionH;
0352       edm::Handle<hgcal::RecoToSimCollectionWithSimClusters> recotosimCollectionH;
0353       event.getByToken(associatorMapRtSim, recotosimCollectionH);
0354       auto recSimColl = *recotosimCollectionH;
0355 
0356       histoProducerAlgo_->fill_simClusterAssociation_histos(histograms.histoProducerAlgo,
0357                                                             ws,
0358                                                             clusterHandle,
0359                                                             clusters,
0360                                                             simClustersHandle,
0361                                                             simClusters,
0362                                                             sCIndices,
0363                                                             inputClusterMask,
0364                                                             *hitMap,
0365                                                             totallayers_to_monitor_,
0366                                                             recSimColl,
0367                                                             simRecColl);
0368 
0369       //General Info on simClusters
0370       LogTrace("HGCalValidator") << "\n# of SimClusters: " << nSimClusters
0371                                  << ", layerClusters mask label: " << label_clustersmask[ws].label() << "\n";
0372     }  //end of loop overs masks
0373   }
0374 
0375   // ##############################################
0376   // Fill layercluster histograms
0377   // ##############################################
0378   int w = 0;  //counter counting the number of sets of histograms
0379   if (doLayerClustersPlots_) {
0380     histoProducerAlgo_->fill_generic_cluster_histos(histograms.histoProducerAlgo,
0381                                                     w,
0382                                                     clusterHandle,
0383                                                     clusters,
0384                                                     caloParticleHandle,
0385                                                     caloParticles,
0386                                                     cPIndices,
0387                                                     selected_cPeff,
0388                                                     *hitMap,
0389                                                     cummatbudg,
0390                                                     totallayers_to_monitor_,
0391                                                     thicknesses_to_monitor_,
0392                                                     recSimColl,
0393                                                     simRecColl);
0394 
0395     for (unsigned int layerclusterIndex = 0; layerclusterIndex < clusters.size(); layerclusterIndex++) {
0396       histoProducerAlgo_->fill_cluster_histos(histograms.histoProducerAlgo, w, clusters[layerclusterIndex]);
0397     }
0398 
0399     //General Info on hgcalLayerClusters
0400     LogTrace("HGCalValidator") << "\n# of layer clusters with " << label_lcl.process() << ":" << label_lcl.label()
0401                                << ":" << label_lcl.instance() << ": " << clusters.size() << "\n";
0402   }
0403 
0404   // ##############################################
0405   // Fill Trackster histograms
0406   // ##############################################
0407   for (unsigned int wml = 0; wml < label_tstTokens.size(); wml++) {
0408     if (doTrackstersPlots_) {
0409       edm::Handle<ticl::TracksterCollection> tracksterHandle;
0410       event.getByToken(label_tstTokens[wml], tracksterHandle);
0411       const ticl::TracksterCollection& tracksters = *tracksterHandle;
0412 
0413       //General Info on Tracksters
0414       LogTrace("HGCalValidator") << "\n# of Tracksters from " << label_tst[wml].process() << ":"
0415                                  << label_tst[wml].label() << ":" << label_tst[wml].instance() << ": "
0416                                  << tracksters.size() << "\n"
0417                                  << std::endl;
0418 
0419       histoProducerAlgo_->fill_trackster_histos(histograms.histoProducerAlgo,
0420                                                 wml,
0421                                                 tracksters,
0422                                                 clusters,
0423                                                 simTracksters,
0424                                                 simTrackstersFromCPs,
0425                                                 cpToSc_SimTrackstersMap,
0426                                                 simClusters,
0427                                                 caloParticleHandle.id(),
0428                                                 caloParticles,
0429                                                 cPIndices,
0430                                                 selected_cPeff,
0431                                                 *hitMap,
0432                                                 totallayers_to_monitor_);
0433     }
0434   }  //end of loop over Trackster input labels
0435 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!