TICLCandidate validation (cms-sw#21)

* changes to simTICLCandidates pdgId

- keep original pdgId for charged simTICLCandidates w/o reco track
- add check for pi0 pdgId

* add ticlCandidates validation

RecoHGCal/TICL/plugins/SimTrackstersProducer.cc

@@ -459,9 +459,7 @@ void SimTrackstersProducer::produce(edm::Event& evt, const edm::EventSetup& es)
 
     auto pdgId = cp.pdgId();
     auto charge = cp.charge();
-    if (cand.trackPtr().isNonnull() and charge == 0) {
-    }
-    if (cand.trackPtr().isNonnull() and charge != 0) {
+    if (cand.trackPtr().isNonnull()) {
       auto const& track = cand.trackPtr().get();
       if (std::abs(pdgId) == 13) {
         cand.setPdgId(pdgId);
@@ -475,7 +473,14 @@ void SimTrackstersProducer::produce(edm::Event& evt, const edm::EventSetup& es)
                                  regressedEnergy);
       cand.setP4(p4);
     } else {  // neutral candidates
-      cand.setPdgId(isHad(pdgId) ? 130 : 22);
+      // a neutral candidate with a charged CaloParticle is charged without a reco track associated with it
+      // set the charge = 0, but keep the real pdgId to keep track of that
+      if (charge != 0)
+        cand.setPdgId(isHad(pdgId) ? 211 : 11);
+      else if (pdgId == 111)
+        cand.setPdgId(pdgId);
+      else
+        cand.setPdgId(isHad(pdgId) ? 130 : 22);
       cand.setCharge(0);
 
       auto particleType = tracksterParticleTypeFromPdgId(cand.pdgId(), 1);

Validation/HGCalValidation/interface/HGCalValidator.h

@@ -22,6 +22,7 @@
 
 #include "DQMServices/Core/interface/DQMGlobalEDAnalyzer.h"
 
+#include "Validation/HGCalValidation/interface/TICLCandidateValidator.h"
 #include "Validation/HGCalValidation/interface/HGVHistoProducerAlgo.h"
 #include "Validation/HGCalValidation/interface/CaloParticleSelector.h"
 #include "RecoLocalCalo/HGCalRecProducers/interface/HGCalClusteringAlgoBase.h"
@@ -75,6 +76,8 @@ class HGCalValidator : public DQMGlobalEDAnalyzer<HGCalValidatorHistograms> {
   const bool doTrackstersPlots_;
   std::string label_TS_, label_TSToCPLinking_, label_TSToSTSPR_;
   std::vector<edm::InputTag> label_clustersmask;
+  const bool doCandidatesPlots_;
+  edm::InputTag label_candidates_;
   const edm::FileInPath cummatbudinxo_;
 
   std::vector<edm::EDGetTokenT<reco::CaloClusterCollection>> labelToken;
@@ -96,6 +99,7 @@ class HGCalValidator : public DQMGlobalEDAnalyzer<HGCalValidatorHistograms> {
   std::unique_ptr<HGVHistoProducerAlgo> histoProducerAlgo_;
 
 private:
+  mutable TICLCandidateValidator candidateVal;
   CaloParticleSelector cpSelector;
   std::shared_ptr<hgcal::RecHitTools> tools_;
   std::map<double, double> cummatbudg;

Validation/HGCalValidation/interface/TICLCandidateValidator.h

@@ -0,0 +1,118 @@
+#ifndef Validation_HGCalValidation_TICLCandidateValidator_h
+#define Validation_HGCalValidation_TICLCandidateValidator_h
+
+#include <iostream>
+#include <vector>
+#include <unordered_map>
+
+#include "FWCore/Framework/interface/ConsumesCollector.h"
+#include "FWCore/ParameterSet/interface/ParameterSet.h"
+#include "FWCore/Framework/interface/EventSetup.h"
+#include "FWCore/Framework/interface/Event.h"
+#include "FWCore/Utilities/interface/EDGetToken.h"
+#include "FWCore/Framework/interface/Frameworkfwd.h"
+
+#include "DataFormats/TrackReco/interface/Track.h"
+#include "DataFormats/HGCalReco/interface/Trackster.h"
+#include "DataFormats/HGCalReco/interface/TICLCandidate.h"
+
+#include "SimDataFormats/Associations/interface/TracksterToSimTracksterHitLCAssociator.h"
+
+#include "DQMServices/Core/interface/DQMStore.h"
+
+class TICLCandidateValidator {
+public:
+  typedef dqm::legacy::DQMStore DQMStore;
+  typedef dqm::legacy::MonitorElement MonitorElement;
+
+  TICLCandidateValidator(){};
+  TICLCandidateValidator(edm::EDGetTokenT<std::vector<TICLCandidate>> TICLCandidates,
+                         edm::EDGetTokenT<std::vector<TICLCandidate>> simTICLCandidatesToken,
+                         edm::EDGetTokenT<std::vector<reco::Track>> recoTracksToken,
+                         edm::EDGetTokenT<std::vector<ticl::Trackster>> trackstersToken,
+                         edm::EDGetTokenT<hgcal::RecoToSimCollectionSimTracksters> associatorMapRtSToken,
+                         edm::EDGetTokenT<hgcal::SimToRecoCollectionSimTracksters> associatorMapStRToken);
+  ~TICLCandidateValidator();
+
+  void bookCandidatesHistos(DQMStore::IBooker& ibook, std::string baseDir);
+
+  void fillCandidateHistos(const edm::Event& event, edm::Handle<ticl::TracksterCollection> simTrackstersCP_h);
+
+private:
+  dqm::reco::MonitorElement* h_tracksters_in_candidate;
+  dqm::reco::MonitorElement* h_candidate_raw_energy;
+  dqm::reco::MonitorElement* h_candidate_regressed_energy;
+  dqm::reco::MonitorElement* h_candidate_pT;
+  dqm::reco::MonitorElement* h_candidate_charge;
+  dqm::reco::MonitorElement* h_candidate_pdgId;
+  dqm::reco::MonitorElement* h_candidate_partType;
+
+  std::vector<dqm::reco::MonitorElement*> h_den_chg_energy_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_energy_candidate_track;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_energy_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_energy_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_chg_pt_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_pt_candidate_track;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_pt_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_pt_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_chg_eta_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_eta_candidate_track;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_eta_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_eta_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_chg_phi_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_phi_candidate_track;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_phi_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_chg_phi_candidate_energy;
+
+  std::vector<dqm::reco::MonitorElement*> h_den_neut_energy_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_neut_energy_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_neut_energy_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_neut_pt_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_neut_pt_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_neut_pt_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_neut_eta_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_neut_eta_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_neut_eta_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_neut_phi_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_neut_phi_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_neut_phi_candidate_energy;
+
+  std::vector<dqm::reco::MonitorElement*> h_den_fake_chg_energy_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_energy_candidate_track;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_energy_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_energy_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_fake_chg_pt_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_pt_candidate_track;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_pt_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_pt_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_fake_chg_eta_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_eta_candidate_track;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_eta_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_eta_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_fake_chg_phi_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_phi_candidate_track;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_phi_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_chg_phi_candidate_energy;
+
+  std::vector<dqm::reco::MonitorElement*> h_den_fake_neut_energy_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_neut_energy_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_neut_energy_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_fake_neut_pt_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_neut_pt_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_neut_pt_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_fake_neut_eta_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_neut_eta_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_neut_eta_candidate_energy;
+  std::vector<dqm::reco::MonitorElement*> h_den_fake_neut_phi_candidate;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_neut_phi_candidate_pdgId;
+  std::vector<dqm::reco::MonitorElement*> h_num_fake_neut_phi_candidate_energy;
+
+  edm::EDGetTokenT<std::vector<TICLCandidate>> TICLCandidatesToken_;
+  edm::EDGetTokenT<std::vector<TICLCandidate>> simTICLCandidatesToken_;
+  edm::EDGetTokenT<std::vector<reco::Track>> recoTracksToken_;
+  edm::EDGetTokenT<std::vector<ticl::Trackster>> trackstersToken_;
+  edm::EDGetTokenT<hgcal::RecoToSimCollectionSimTracksters> associatorMapRtSToken_;
+  edm::EDGetTokenT<hgcal::SimToRecoCollectionSimTracksters> associatorMapStRToken_;
+};
+
+#endif

Validation/HGCalValidation/plugins/HGCalValidator.cc

@@ -30,6 +30,8 @@ HGCalValidator::HGCalValidator(const edm::ParameterSet& pset)
       label_TSToCPLinking_(pset.getParameter<std::string>("label_TSToCPLinking")),
       label_TSToSTSPR_(pset.getParameter<std::string>("label_TSToSTSPR")),
       label_clustersmask(pset.getParameter<std::vector<edm::InputTag>>("LayerClustersInputMask")),
+      doCandidatesPlots_(pset.getUntrackedParameter<bool>("doCandidatesPlots")),
+      label_candidates_(pset.getParameter<edm::InputTag>("ticlCandidates")),
       cummatbudinxo_(pset.getParameter<edm::FileInPath>("cummatbudinxo")) {
   //In this way we can easily generalize to associations between other objects also.
   const edm::InputTag& label_cp_effic_tag = pset.getParameter<edm::InputTag>("label_cp_effic");
@@ -55,6 +57,29 @@ HGCalValidator::HGCalValidator(const edm::ParameterSet& pset)
 
   layerclusters_ = consumes<reco::CaloClusterCollection>(label_lcl);
 
+  if (doCandidatesPlots_) {
+    edm::EDGetTokenT<std::vector<TICLCandidate>> TICLCandidatesToken =
+        consumes<std::vector<TICLCandidate>>(pset.getParameter<edm::InputTag>("ticlTrackstersMerge"));
+    edm::EDGetTokenT<std::vector<TICLCandidate>> simTICLCandidatesToken =
+        consumes<std::vector<TICLCandidate>>(pset.getParameter<edm::InputTag>("simTiclCandidates"));
+    edm::EDGetTokenT<std::vector<reco::Track>> recoTracksToken =
+        consumes<std::vector<reco::Track>>(pset.getParameter<edm::InputTag>("recoTracks"));
+    edm::EDGetTokenT<std::vector<ticl::Trackster>> trackstersToken =
+        consumes<std::vector<ticl::Trackster>>(pset.getParameter<edm::InputTag>("ticlTrackstersMerge"));
+    //consumes<std::vector<ticl::Trackster>>(pset.getParameter<edm::InputTag>("trackstersclue3d"));
+    edm::EDGetTokenT<hgcal::RecoToSimCollectionSimTracksters> associatorMapRtSToken =
+        consumes<hgcal::SimToRecoCollectionSimTracksters>(pset.getParameter<edm::InputTag>("mergeRecoToSimAssociator"));
+    edm::EDGetTokenT<hgcal::SimToRecoCollectionSimTracksters> associatorMapStRToken =
+        consumes<hgcal::SimToRecoCollectionSimTracksters>(pset.getParameter<edm::InputTag>("mergeSimToRecoAssociator"));
+
+    candidateVal = TICLCandidateValidator(TICLCandidatesToken,
+                                          simTICLCandidatesToken,
+                                          recoTracksToken,
+                                          trackstersToken,
+                                          associatorMapRtSToken,
+                                          associatorMapStRToken);
+  }
+
   for (auto& itag : label_tst) {
     label_tstTokens.push_back(consumes<ticl::TracksterCollection>(itag));
   }
@@ -219,6 +244,13 @@ void HGCalValidator::bookHistograms(DQMStore::IBooker& ibook,
           ibook, histograms.histoProducerAlgo, HGVHistoProducerAlgo::validationType::PatternRecognition);
     }
   }  //end of booking Tracksters loop
+
+  // Booking histograms concerning TICL candidates
+  if (doCandidatesPlots_) {
+    ibook.cd();
+    ibook.setCurrentFolder(dirName_ + label_candidates_.label());
+    candidateVal.bookCandidatesHistos(ibook, dirName_ + label_candidates_.label());
+  }
 }
 
 void HGCalValidator::cpParametersAndSelection(const Histograms& histograms,
@@ -432,4 +464,9 @@ void HGCalValidator::dqmAnalyze(const edm::Event& event,
                                                 totallayers_to_monitor_);
     }
   }  //end of loop over Trackster input labels
+
+  // tracksters histograms
+  if (doCandidatesPlots_) {
+    candidateVal.fillCandidateHistos(event, simTracksterFromCPHandle);
+  }
 }

Validation/HGCalValidation/python/HGCalPostProcessor_cff.py

@@ -3,7 +3,7 @@
 from Validation.HGCalValidation.HGCalSimHitsClient_cff import *
 from Validation.HGCalValidation.HGCalDigiClient_cff    import *
 from Validation.HGCalValidation.HGCalRecHitsClient_cff import *
-from Validation.HGCalValidation.PostProcessorHGCAL_cfi import postProcessorHGCALlayerclusters,postProcessorHGCALsimclusters,postProcessorHGCALTracksters
+from Validation.HGCalValidation.PostProcessorHGCAL_cfi import postProcessorHGCALlayerclusters,postProcessorHGCALsimclusters,postProcessorHGCALTracksters,postProcessorHGCALCandidates
 
 hgcalPostProcessor = cms.Sequence(hgcalSimHitClientEE
     + hgcalSimHitClientHEF
@@ -18,4 +18,5 @@
 hgcalValidatorPostProcessor = cms.Sequence(
     postProcessorHGCALlayerclusters+
     postProcessorHGCALsimclusters+
-    postProcessorHGCALTracksters)
+    postProcessorHGCALTracksters+
+    postProcessorHGCALCandidates)

Validation/HGCalValidation/python/HGCalValidator_cfi.py

@@ -51,6 +51,16 @@
     label_TS = cms.string("Morphology"),
     label_TSToCPLinking = cms.string("TSToCP_linking"),
     label_TSToSTSPR = cms.string("TSToSTS_patternRecognition"),
+    #candidates plots
+    doCandidatesPlots = cms.untracked.bool(True),
+    ticlCandidates = cms.string("ticlCandidates"),
+
+    ticlTrackstersMerge = cms.InputTag("ticlTrackstersMerge"),
+    simTiclCandidates = cms.InputTag("ticlSimTracksters"),
+    recoTracks = cms.InputTag("generalTracks"),
+    trackstersclue3d = cms.InputTag("ticlTrackstersCLUE3DHigh"),
+    mergeRecoToSimAssociator = cms.InputTag("tracksterSimTracksterAssociationLinking", "recoToSim"),
+    mergeSimToRecoAssociator = cms.InputTag("tracksterSimTracksterAssociationLinking", "simToReco"),
 
     #The cumulative material budget in front of each layer. To be more specific, it
     #is the material budget just in front of the active material (not including it).

Validation/HGCalValidation/python/PostProcessorHGCAL_cfi.py

@@ -82,3 +82,39 @@
   outputFileName = cms.untracked.string(""),
   verbose = cms.untracked.uint32(4)
 )
+
+neutrals = ["photons", "neutral_pions", "neutral_hadrons"]
+charged = ["electrons", "muons", "charged_hadrons"]
+variables = ["energy", "pt", "eta", "phi"]
+subDirsCandidates = [prefix + hgcalValidator.ticlCandidates.value() + "/" + c for cands in (neutrals, charged) for c in cands]
+eff_candidates = []
+
+for c in charged:
+    for var in variables:
+        eff_candidates.append("eff_"+c+"_track_"+var+" '"+c.replace("_", " ")+" candidates track efficiency vs "+var+"' num_track_cand_vs_"+var+"_"+c+" den_cand_vs_"+var+"_"+c)
+        eff_candidates.append("eff_"+c+"_pid_"+var+" '"+c.replace("_", " ")+" candidates track + pid efficiency vs "+var+"' num_pid_cand_vs_"+var+"_"+c+" den_cand_vs_"+var+"_"+c)
+        eff_candidates.append("eff_"+c+"_energy_"+var+" '"+c.replace("_", " ")+" candidates track + pid + energy efficiency vs "+var+"' num_energy_cand_vs_"+var+"_"+c+" den_cand_vs_"+var+"_"+c)
+for n in neutrals:
+    for var in variables:
+        eff_candidates.append("eff_"+n+"_pid_"+var+" '"+n.replace("_", " ")+" candidates pid efficiency vs "+var+"' num_pid_cand_vs_"+var+"_"+n+" den_cand_vs_"+var+"_"+n)
+        eff_candidates.append("eff_"+n+"_energy_"+var+" '"+n.replace("_", " ")+" candidates pid + energy efficiency vs "+var+"' num_energy_cand_vs_"+var+"_"+n+" den_cand_vs_"+var+"_"+n)
+
+for c in charged:
+    for var in variables:
+        eff_candidates.append("fake_"+c+"_track_"+var+" '"+c.replace("_", " ")+" candidates track fake vs "+var+"' num_fake_track_cand_vs_"+var+"_"+c+" den_fake_cand_vs_"+var+"_"+c)
+        eff_candidates.append("fake_"+c+"_pid_"+var+" '"+c.replace("_", " ")+" candidates track + pid fake vs "+var+"' num_fake_pid_cand_vs_"+var+"_"+c+" den_fake_cand_vs_"+var+"_"+c)
+        eff_candidates.append("fake_"+c+"_energy_"+var+" '"+c.replace("_", " ")+" candidates track + pid + energy fake vs "+var+"' num_fake_energy_cand_vs_"+var+"_"+c+" den_fake_cand_vs_"+var+"_"+c)
+for n in neutrals:
+    for var in variables:
+        eff_candidates.append("fake_"+n+"_pid_"+var+" '"+n.replace("_", " ")+" candidates pid fake vs "+var+"' num_fake_pid_cand_vs_"+var+"_"+n+" den_fake_cand_vs_"+var+"_"+n)
+        eff_candidates.append("fake_"+n+"_energy_"+var+" '"+n.replace("_", " ")+" candidates pid + energy fake vs "+var+"' num_fake_energy_cand_vs_"+var+"_"+n+" den_fake_cand_vs_"+var+"_"+n)
+
+postProcessorHGCALCandidates = DQMEDHarvester('DQMGenericClient',
+  subDirs = cms.untracked.vstring(subDirsCandidates),
+  efficiency = cms.vstring(eff_candidates),
+  resolution = cms.vstring(),
+  cumulativeDists = cms.untracked.vstring(),
+  noFlowDists = cms.untracked.vstring(),
+  outputFileName = cms.untracked.string(""),
+  verbose = cms.untracked.uint32(4)
+)