CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
GsfElectronProducer.cc
Go to the documentation of this file.
27 
28 #include <array>
29 
30 using namespace reco;
31 
32 namespace {
33 
34  void setMVAOutputs(reco::GsfElectronCollection& electrons,
37  bool dnnPFidEnabled,
38  const std::vector<tensorflow::Session*>& tfSessions) {
39  std::vector<GsfElectron::MvaOutput> mva_outputs(electrons.size());
40  size_t iele = 0;
41  for (auto& el : electrons) {
42  GsfElectron::MvaOutput mvaOutput;
43  mvaOutput.mva_e_pi = hoc->sElectronMVAEstimator->mva(el, vertices);
44  mvaOutput.mva_Isolated = hoc->iElectronMVAEstimator->mva(el, vertices.size());
45  if (dnnPFidEnabled) {
46  mva_outputs[iele] = mvaOutput;
47  } else {
48  el.setMvaOutput(mvaOutput);
49  }
50  iele++;
51  }
52  if (dnnPFidEnabled) {
53  // Here send the list of electrons to the ElectronDNNEstimator and get back the values for all the electrons in one go
54  LogDebug("GsfElectronProducer") << "Getting DNN PFId for ele";
55  const auto& dnn_ele_pfid = hoc->iElectronDNNEstimator->evaluate(electrons, tfSessions);
56  int jele = 0;
57  for (auto& el : electrons) {
58  const auto& values = dnn_ele_pfid[jele];
59  // get the previous values
60  auto& mvaOutput = mva_outputs[jele];
61  mvaOutput.dnn_e_sigIsolated = values[0];
62  mvaOutput.dnn_e_sigNonIsolated = values[1];
63  mvaOutput.dnn_e_bkgNonIsolated = values[2];
64  mvaOutput.dnn_e_bkgTau = values[3];
65  mvaOutput.dnn_e_bkgPhoton = values[4];
66  el.setMvaOutput(mvaOutput);
67  jele++;
68  }
69  }
70  }
71 
72  // Something more clever has to be found. The collections are small, so the timing is not
73  // an issue here; but it is clearly suboptimal
74 
75  auto matchWithPFCandidates(std::vector<reco::PFCandidate> const& pfCandidates) {
76  std::map<reco::GsfTrackRef, reco::GsfElectron::MvaInput> gsfMVAInputs{};
77 
78  //Loop over the collection of PFFCandidates
79  for (auto const& pfCand : pfCandidates) {
81  // First check that the GsfTrack is non null
82  if (pfCand.gsfTrackRef().isNonnull()) {
84  input.earlyBrem = pfCand.egammaExtraRef()->mvaVariable(reco::PFCandidateEGammaExtra::MVA_FirstBrem);
85  input.lateBrem = pfCand.egammaExtraRef()->mvaVariable(reco::PFCandidateEGammaExtra::MVA_LateBrem);
86  input.deltaEta = pfCand.egammaExtraRef()->mvaVariable(reco::PFCandidateEGammaExtra::MVA_DeltaEtaTrackCluster);
87  input.sigmaEtaEta = pfCand.egammaExtraRef()->sigmaEtaEta();
88  input.hadEnergy = pfCand.egammaExtraRef()->hadEnergy();
89  gsfMVAInputs[pfCand.gsfTrackRef()] = input;
90  }
91  }
92  return gsfMVAInputs;
93  }
94 
95  void logElectrons(reco::GsfElectronCollection const& electrons, edm::Event const& event, const std::string& title) {
96  LogTrace("GsfElectronAlgo") << "========== " << title << " ==========";
97  LogTrace("GsfElectronAlgo") << "Event: " << event.id();
98  LogTrace("GsfElectronAlgo") << "Number of electrons: " << electrons.size();
99  for (auto const& ele : electrons) {
100  LogTrace("GsfElectronAlgo") << "Electron with charge, pt, eta, phi: " << ele.charge() << " , " << ele.pt()
101  << " , " << ele.eta() << " , " << ele.phi();
102  }
103  LogTrace("GsfElectronAlgo") << "=================================================";
104  }
105 
106 } // namespace
107 
108 class GsfElectronProducer : public edm::stream::EDProducer<edm::GlobalCache<GsfElectronAlgo::HeavyObjectCache>> {
109 public:
111 
113 
114  static std::unique_ptr<GsfElectronAlgo::HeavyObjectCache> initializeGlobalCache(const edm::ParameterSet& conf) {
115  return std::make_unique<GsfElectronAlgo::HeavyObjectCache>(conf);
116  }
117 
118  void endStream() override;
119 
121 
122  // ------------ method called to produce the data ------------
123  void produce(edm::Event& event, const edm::EventSetup& setup) override;
124 
125 private:
126  std::unique_ptr<GsfElectronAlgo> algo_;
127 
128  // configurables
133 
135 
136  bool isPreselected(reco::GsfElectron const& ele) const;
137  void setAmbiguityData(reco::GsfElectronCollection& electrons,
138  edm::Event const& event,
139  bool ignoreNotPreselected = true) const;
140 
141  // check expected configuration of previous modules
143  void checkEcalSeedingParameters(edm::ParameterSet const&);
144 
148 
149  const bool useGsfPfRecTracks_;
150 
152 
154 
155  std::vector<tensorflow::Session*> tfSessions_;
156 };
157 
160  // input collections
161  desc.add<edm::InputTag>("gsfElectronCoresTag", {"ecalDrivenGsfElectronCores"});
162  desc.add<edm::InputTag>("vtxTag", {"offlinePrimaryVertices"});
163  desc.add<edm::InputTag>("conversionsTag", {"allConversions"});
164  desc.add<edm::InputTag>("gsfPfRecTracksTag", {"pfTrackElec"});
165  desc.add<edm::InputTag>("barrelRecHitCollectionTag", {"ecalRecHit", "EcalRecHitsEB"});
166  desc.add<edm::InputTag>("endcapRecHitCollectionTag", {"ecalRecHit", "EcalRecHitsEE"});
167  desc.add<edm::InputTag>("seedsTag", {"ecalDrivenElectronSeeds"});
168  desc.add<edm::InputTag>("beamSpotTag", {"offlineBeamSpot"});
169  desc.add<edm::InputTag>("egmPFCandidatesTag", {"particleFlowEGamma"});
170 
171  // steering
172  desc.add<bool>("useDefaultEnergyCorrection", true);
173  desc.add<bool>("useCombinationRegression", false);
174  desc.add<bool>("ecalDrivenEcalEnergyFromClassBasedParameterization", true);
175  desc.add<bool>("ecalDrivenEcalErrorFromClassBasedParameterization", true);
176  desc.add<bool>("applyPreselection", false);
177  desc.add<bool>("useEcalRegression", false);
178  desc.add<bool>("applyAmbResolution", false);
179  desc.add<bool>("ignoreNotPreselected", true);
180  desc.add<bool>("useGsfPfRecTracks", true);
181  desc.add<bool>("pureTrackerDrivenEcalErrorFromSimpleParameterization", true);
182  desc.add<unsigned int>("ambSortingStrategy", 1);
183  desc.add<unsigned int>("ambClustersOverlapStrategy", 1);
184  desc.add<bool>("fillConvVtxFitProb", true);
185  desc.add<bool>("resetMvaValuesUsingPFCandidates", false);
186 
187  // Ecal rec hits configuration
188  desc.add<std::vector<std::string>>("recHitFlagsToBeExcludedBarrel");
189  desc.add<std::vector<std::string>>("recHitFlagsToBeExcludedEndcaps");
190  desc.add<std::vector<std::string>>("recHitSeverityToBeExcludedBarrel");
191  desc.add<std::vector<std::string>>("recHitSeverityToBeExcludedEndcaps");
192 
193  // Hcal rec hits configuration
194  desc.add<bool>("checkHcalStatus", true);
195  desc.add<edm::InputTag>("hbheRecHits", edm::InputTag("hbhereco"));
196  desc.add<std::vector<double>>("recHitEThresholdHB", {0., 0., 0., 0.});
197  desc.add<std::vector<double>>("recHitEThresholdHE", {0., 0., 0., 0., 0., 0., 0.});
198  desc.add<int>("maxHcalRecHitSeverity", 999999);
199  desc.add<bool>("hcalRun2EffDepth", false);
200 
201  // Isolation algos configuration
202  desc.add("trkIsol03Cfg", EleTkIsolFromCands::pSetDescript());
203  desc.add("trkIsol04Cfg", EleTkIsolFromCands::pSetDescript());
204  desc.add("trkIsolHEEP03Cfg", EleTkIsolFromCands::pSetDescript());
205  desc.add("trkIsolHEEP04Cfg", EleTkIsolFromCands::pSetDescript());
206  desc.add<bool>("useNumCrystals", true);
207  desc.add<double>("etMinBarrel", 0.0);
208  desc.add<double>("etMinEndcaps", 0.11);
209  desc.add<double>("etMinHcal", 0.0);
210  desc.add<double>("eMinBarrel", 0.095);
211  desc.add<double>("eMinEndcaps", 0.0);
212  desc.add<double>("intRadiusEcalBarrel", 3.0);
213  desc.add<double>("intRadiusEcalEndcaps", 3.0);
214  desc.add<double>("intRadiusHcal", 0.15);
215  desc.add<double>("jurassicWidth", 1.5);
216  desc.add<bool>("vetoClustered", false);
217 
218  // backward compatibility mechanism for ctf tracks
219  desc.add<bool>("ctfTracksCheck", true);
220  desc.add<edm::InputTag>("ctfTracksTag", {"generalTracks"});
221 
222  desc.add<double>("MaxElePtForOnlyMVA", 50.0);
223  desc.add<double>("PreSelectMVA", -0.1);
224 
225  {
227  psd0.add<double>("minSCEtBarrel", 4.0);
228  psd0.add<double>("minSCEtEndcaps", 4.0);
229  psd0.add<double>("minEOverPBarrel", 0.0);
230  psd0.add<double>("minEOverPEndcaps", 0.0);
231  psd0.add<double>("maxEOverPBarrel", 999999999.0);
232  psd0.add<double>("maxEOverPEndcaps", 999999999.0);
233  psd0.add<double>("maxDeltaEtaBarrel", 0.02);
234  psd0.add<double>("maxDeltaEtaEndcaps", 0.02);
235  psd0.add<double>("maxDeltaPhiBarrel", 0.15);
236  psd0.add<double>("maxDeltaPhiEndcaps", 0.15);
237  psd0.add<double>("hOverEConeSize", 0.15);
238  psd0.add<double>("maxHOverEBarrelCone", 0.15);
239  psd0.add<double>("maxHOverEEndcapsCone", 0.15);
240  psd0.add<double>("maxHBarrelCone", 0.0);
241  psd0.add<double>("maxHEndcapsCone", 0.0);
242  psd0.add<double>("maxHOverEBarrelBc", 0.15);
243  psd0.add<double>("maxHOverEEndcapsBc", 0.15);
244  psd0.add<double>("maxHBarrelBc", 0.0);
245  psd0.add<double>("maxHEndcapsBc", 0.0);
246  psd0.add<double>("maxSigmaIetaIetaBarrel", 999999999.0);
247  psd0.add<double>("maxSigmaIetaIetaEndcaps", 999999999.0);
248  psd0.add<double>("maxFbremBarrel", 999999999.0);
249  psd0.add<double>("maxFbremEndcaps", 999999999.0);
250  psd0.add<bool>("isBarrel", false);
251  psd0.add<bool>("isEndcaps", false);
252  psd0.add<bool>("isFiducial", false);
253  psd0.add<bool>("seedFromTEC", true);
254  psd0.add<double>("maxTIP", 999999999.0);
255  psd0.add<double>("multThresEB", EgammaLocalCovParamDefaults::kMultThresEB);
256  psd0.add<double>("multThresEE", EgammaLocalCovParamDefaults::kMultThresEE);
257  // preselection parameters
258  desc.add<edm::ParameterSetDescription>("preselection", psd0);
259  }
260 
261  // Corrections
262  desc.add<std::string>("crackCorrectionFunction", "EcalClusterCrackCorrection");
263 
264  desc.add<bool>("ecalWeightsFromDB", true);
265  desc.add<std::vector<std::string>>("ecalRefinedRegressionWeightFiles", {})
266  ->setComment("if not from DB. Otherwise, keep empty");
267  desc.add<bool>("combinationWeightsFromDB", true);
268  desc.add<std::vector<std::string>>("combinationRegressionWeightFile", {})
269  ->setComment("if not from DB. Otherwise, keep empty");
270 
271  // regression. The labels are needed in all cases.
272  desc.add<std::vector<std::string>>("ecalRefinedRegressionWeightLabels", {});
273  desc.add<std::vector<std::string>>("combinationRegressionWeightLabels", {});
274 
275  desc.add<std::vector<std::string>>(
276  "ElecMVAFilesString",
277  {
278  "RecoEgamma/ElectronIdentification/data/TMVA_Category_BDTSimpleCat_10_17Feb2011.weights.xml",
279  "RecoEgamma/ElectronIdentification/data/TMVA_Category_BDTSimpleCat_12_17Feb2011.weights.xml",
280  "RecoEgamma/ElectronIdentification/data/TMVA_Category_BDTSimpleCat_20_17Feb2011.weights.xml",
281  "RecoEgamma/ElectronIdentification/data/TMVA_Category_BDTSimpleCat_22_17Feb2011.weights.xml",
282  });
283  desc.add<std::vector<std::string>>(
284  "SoftElecMVAFilesString",
285  {
286  "RecoEgamma/ElectronIdentification/data/TMVA_BDTSoftElectrons_7Feb2014.weights.xml",
287  });
288 
289  {
291  psd1.add<bool>("enabled", false);
292  psd1.add<std::string>("inputTensorName", "FirstLayer_input");
293  psd1.add<std::string>("outputTensorName", "sequential/FinalLayer/Softmax");
294  psd1.add<uint>("outputDim", 5); // Number of output nodes of DNN
295  psd1.add<std::vector<std::string>>(
296  "modelsFiles",
297  {"RecoEgamma/ElectronIdentification/data/Ele_PFID_dnn/v1/lowpT/lowpT_modelDNN.pb",
298  "RecoEgamma/ElectronIdentification/data/Ele_PFID_dnn/v1/highpTEB/highpTEB_modelDNN.pb",
299  "RecoEgamma/ElectronIdentification/data/Ele_PFID_dnn/v1/highpTEE/highpTEE_modelDNN.pb"});
300  psd1.add<std::vector<std::string>>(
301  "scalersFiles",
302  {"RecoEgamma/ElectronIdentification/data/Ele_PFID_dnn/v1/lowpT/lowpT_scaler.txt",
303  "RecoEgamma/ElectronIdentification/data/Ele_PFID_dnn/v1/highpTEB/highpTEB_scaler.txt",
304  "RecoEgamma/ElectronIdentification/data/Ele_PFID_dnn/v1/highpTEE/highpTEE_scaler.txt"});
305  psd1.add<bool>("useEBModelInGap", true);
306  // preselection parameters
307  desc.add<edm::ParameterSetDescription>("EleDNNPFid", psd1);
308  }
309 
312  {
314  psd0.add<edm::InputTag>("pfClusterProducer", edm::InputTag("particleFlowClusterECAL"));
315  psd0.add<double>("drMax", 0.3);
316  psd0.add<double>("drVetoBarrel", 0.0);
317  psd0.add<double>("drVetoEndcap", 0.0);
318  psd0.add<double>("etaStripBarrel", 0.0);
319  psd0.add<double>("etaStripEndcap", 0.0);
320  psd0.add<double>("energyBarrel", 0.0);
321  psd0.add<double>("energyEndcap", 0.0);
322  desc.add<edm::ParameterSetDescription>("pfECALClusIsolCfg", psd0);
323  }
324 
326  {
328  psd0.add<edm::InputTag>("pfClusterProducerHCAL", edm::InputTag("particleFlowClusterHCAL"));
329  psd0.add<edm::InputTag>("pfClusterProducerHFEM", edm::InputTag(""));
330  psd0.add<edm::InputTag>("pfClusterProducerHFHAD", edm::InputTag(""));
331  psd0.add<bool>("useHF", false);
332  psd0.add<double>("drMax", 0.3);
333  psd0.add<double>("drVetoBarrel", 0.0);
334  psd0.add<double>("drVetoEndcap", 0.0);
335  psd0.add<double>("etaStripBarrel", 0.0);
336  psd0.add<double>("etaStripEndcap", 0.0);
337  psd0.add<double>("energyBarrel", 0.0);
338  psd0.add<double>("energyEndcap", 0.0);
339  psd0.add<bool>("useEt", true);
340  desc.add<edm::ParameterSetDescription>("pfHCALClusIsolCfg", psd0);
341  }
342 
343  descriptions.add("gsfElectronProducerDefault", desc);
344 }
345 
346 namespace {
347  GsfElectronAlgo::CutsConfiguration makeCutsConfiguration(edm::ParameterSet const& pset) {
349  .minSCEtBarrel = pset.getParameter<double>("minSCEtBarrel"),
350  .minSCEtEndcaps = pset.getParameter<double>("minSCEtEndcaps"),
351  .maxEOverPBarrel = pset.getParameter<double>("maxEOverPBarrel"),
352  .maxEOverPEndcaps = pset.getParameter<double>("maxEOverPEndcaps"),
353  .minEOverPBarrel = pset.getParameter<double>("minEOverPBarrel"),
354  .minEOverPEndcaps = pset.getParameter<double>("minEOverPEndcaps"),
355  .maxHOverEBarrelCone = pset.getParameter<double>("maxHOverEBarrelCone"),
356  .maxHOverEEndcapsCone = pset.getParameter<double>("maxHOverEEndcapsCone"),
357  .maxHBarrelCone = pset.getParameter<double>("maxHBarrelCone"),
358  .maxHEndcapsCone = pset.getParameter<double>("maxHEndcapsCone"),
359  .maxHOverEBarrelBc = pset.getParameter<double>("maxHOverEBarrelBc"),
360  .maxHOverEEndcapsBc = pset.getParameter<double>("maxHOverEEndcapsBc"),
361  .maxHBarrelBc = pset.getParameter<double>("maxHBarrelBc"),
362  .maxHEndcapsBc = pset.getParameter<double>("maxHEndcapsBc"),
363  .maxDeltaEtaBarrel = pset.getParameter<double>("maxDeltaEtaBarrel"),
364  .maxDeltaEtaEndcaps = pset.getParameter<double>("maxDeltaEtaEndcaps"),
365  .maxDeltaPhiBarrel = pset.getParameter<double>("maxDeltaPhiBarrel"),
366  .maxDeltaPhiEndcaps = pset.getParameter<double>("maxDeltaPhiEndcaps"),
367  .maxSigmaIetaIetaBarrel = pset.getParameter<double>("maxSigmaIetaIetaBarrel"),
368  .maxSigmaIetaIetaEndcaps = pset.getParameter<double>("maxSigmaIetaIetaEndcaps"),
369  .maxFbremBarrel = pset.getParameter<double>("maxFbremBarrel"),
370  .maxFbremEndcaps = pset.getParameter<double>("maxFbremEndcaps"),
371  .isBarrel = pset.getParameter<bool>("isBarrel"),
372  .isEndcaps = pset.getParameter<bool>("isEndcaps"),
373  .isFiducial = pset.getParameter<bool>("isFiducial"),
374  .maxTIP = pset.getParameter<double>("maxTIP"),
375  .seedFromTEC = pset.getParameter<bool>("seedFromTEC"),
376  .multThresEB = pset.getParameter<double>("multThresEB"),
377  .multThresEE = pset.getParameter<double>("multThresEE"),
378  };
379  }
380 }; // namespace
381 
383  : cutsCfg_{makeCutsConfiguration(cfg.getParameter<edm::ParameterSet>("preselection"))},
384  ecalSeedingParametersChecked_(false),
385  electronPutToken_(produces<GsfElectronCollection>()),
386  gsfPfRecTracksTag_(consumes(cfg.getParameter<edm::InputTag>("gsfPfRecTracksTag"))),
387  useGsfPfRecTracks_(cfg.getParameter<bool>("useGsfPfRecTracks")),
388  resetMvaValuesUsingPFCandidates_(cfg.getParameter<bool>("resetMvaValuesUsingPFCandidates")) {
389  if (resetMvaValuesUsingPFCandidates_) {
390  egmPFCandidateCollection_ = consumes(cfg.getParameter<edm::InputTag>("egmPFCandidatesTag"));
391  }
392 
393  inputCfg_.gsfElectronCores = consumes(cfg.getParameter<edm::InputTag>("gsfElectronCoresTag"));
394  inputCfg_.hbheRecHitsTag = consumes(cfg.getParameter<edm::InputTag>("hbheRecHits"));
395  inputCfg_.barrelRecHitCollection = consumes(cfg.getParameter<edm::InputTag>("barrelRecHitCollectionTag"));
396  inputCfg_.endcapRecHitCollection = consumes(cfg.getParameter<edm::InputTag>("endcapRecHitCollectionTag"));
397  inputCfg_.ctfTracks = consumes(cfg.getParameter<edm::InputTag>("ctfTracksTag"));
398  // used to check config consistency with seeding
399  inputCfg_.seedsTag = consumes(cfg.getParameter<edm::InputTag>("seedsTag"));
400  inputCfg_.beamSpotTag = consumes(cfg.getParameter<edm::InputTag>("beamSpotTag"));
401  inputCfg_.vtxCollectionTag = consumes(cfg.getParameter<edm::InputTag>("vtxTag"));
402  if (cfg.getParameter<bool>("fillConvVtxFitProb"))
403  inputCfg_.conversions = consumes(cfg.getParameter<edm::InputTag>("conversionsTag"));
404 
405  // inputs for PFCluster isolation
406  const edm::ParameterSet& pfECALClusIsolCfg = cfg.getParameter<edm::ParameterSet>("pfECALClusIsolCfg");
407  const edm::ParameterSet& pfHCALClusIsolCfg = cfg.getParameter<edm::ParameterSet>("pfHCALClusIsolCfg");
408  inputCfg_.pfClusterProducer =
409  consumes<reco::PFClusterCollection>(pfECALClusIsolCfg.getParameter<edm::InputTag>("pfClusterProducer"));
410  inputCfg_.pfClusterProducerHCAL = consumes(pfHCALClusIsolCfg.getParameter<edm::InputTag>("pfClusterProducerHCAL"));
411  inputCfg_.pfClusterProducerHFEM = consumes(pfHCALClusIsolCfg.getParameter<edm::InputTag>("pfClusterProducerHFEM"));
412  inputCfg_.pfClusterProducerHFHAD = consumes(pfHCALClusIsolCfg.getParameter<edm::InputTag>("pfClusterProducerHFHAD"));
413 
414  // Config for PFID dnn
415  const auto& pset_dnn = cfg.getParameter<edm::ParameterSet>("EleDNNPFid");
416  dnnPFidEnabled_ = pset_dnn.getParameter<bool>("enabled");
417 
418  strategyCfg_.useDefaultEnergyCorrection = cfg.getParameter<bool>("useDefaultEnergyCorrection");
419 
420  strategyCfg_.applyPreselection = cfg.getParameter<bool>("applyPreselection");
421  strategyCfg_.ecalDrivenEcalEnergyFromClassBasedParameterization =
422  cfg.getParameter<bool>("ecalDrivenEcalEnergyFromClassBasedParameterization");
423  strategyCfg_.ecalDrivenEcalErrorFromClassBasedParameterization =
424  cfg.getParameter<bool>("ecalDrivenEcalErrorFromClassBasedParameterization");
425  strategyCfg_.pureTrackerDrivenEcalErrorFromSimpleParameterization =
426  cfg.getParameter<bool>("pureTrackerDrivenEcalErrorFromSimpleParameterization");
427  strategyCfg_.applyAmbResolution = cfg.getParameter<bool>("applyAmbResolution");
428  strategyCfg_.ignoreNotPreselected = cfg.getParameter<bool>("ignoreNotPreselected");
429  strategyCfg_.ambSortingStrategy = cfg.getParameter<unsigned>("ambSortingStrategy");
430  strategyCfg_.ambClustersOverlapStrategy = cfg.getParameter<unsigned>("ambClustersOverlapStrategy");
431  strategyCfg_.ctfTracksCheck = cfg.getParameter<bool>("ctfTracksCheck");
432  strategyCfg_.PreSelectMVA = cfg.getParameter<double>("PreSelectMVA");
433  strategyCfg_.MaxElePtForOnlyMVA = cfg.getParameter<double>("MaxElePtForOnlyMVA");
434  strategyCfg_.useEcalRegression = cfg.getParameter<bool>("useEcalRegression");
435  strategyCfg_.useCombinationRegression = cfg.getParameter<bool>("useCombinationRegression");
436  strategyCfg_.fillConvVtxFitProb = cfg.getParameter<bool>("fillConvVtxFitProb");
437  strategyCfg_.computePfClusterIso = dnnPFidEnabled_;
438 
439  // hcal helpers
440  auto const& psetPreselection = cfg.getParameter<edm::ParameterSet>("preselection");
441  hcalCfg_.hOverEConeSize = psetPreselection.getParameter<double>("hOverEConeSize");
442  if (hcalCfg_.hOverEConeSize > 0) {
443  hcalCfg_.onlyBehindCluster = false;
444  hcalCfg_.checkHcalStatus = cfg.getParameter<bool>("checkHcalStatus");
445 
446  //hcalCfg_.hbheRecHits = consumes<HBHERecHitCollection>(cfg.getParameter<edm::InputTag>("hbheRecHits"));
447  hcalCfg_.hbheRecHits = consumes<HBHERecHitCollection>(cfg.getParameter<edm::InputTag>("hbheRecHits"));
448 
449  hcalCfg_.eThresHB = cfg.getParameter<EgammaHcalIsolation::arrayHB>("recHitEThresholdHB");
450  hcalCfg_.maxSeverityHB = cfg.getParameter<int>("maxHcalRecHitSeverity");
451  hcalCfg_.eThresHE = cfg.getParameter<EgammaHcalIsolation::arrayHE>("recHitEThresholdHE");
452  hcalCfg_.maxSeverityHE = hcalCfg_.maxSeverityHB;
453  }
454 
455  hcalCfgBc_.hOverEConeSize = 0.;
456  hcalCfgBc_.onlyBehindCluster = true;
457  hcalCfgBc_.checkHcalStatus = cfg.getParameter<bool>("checkHcalStatus");
458 
459  //hcalCfgBc_.hbheRecHits = consumes<HBHERecHitCollection>(cfg.getParameter<edm::InputTag>("hbheRecHits"));
460  hcalCfgBc_.hbheRecHits = consumes<HBHERecHitCollection>(cfg.getParameter<edm::InputTag>("hbheRecHits"));
461 
462  hcalCfgBc_.eThresHB = cfg.getParameter<EgammaHcalIsolation::arrayHB>("recHitEThresholdHB");
463  hcalCfgBc_.maxSeverityHB = cfg.getParameter<int>("maxHcalRecHitSeverity");
464  hcalCfgBc_.eThresHE = cfg.getParameter<EgammaHcalIsolation::arrayHE>("recHitEThresholdHE");
465  hcalCfgBc_.maxSeverityHE = hcalCfgBc_.maxSeverityHB;
466 
467  hcalRun2EffDepth_ = cfg.getParameter<bool>("hcalRun2EffDepth");
468 
469  // Ecal rec hits configuration
471  auto const& flagnamesbarrel = cfg.getParameter<std::vector<std::string>>("recHitFlagsToBeExcludedBarrel");
472  recHitsCfg.recHitFlagsToBeExcludedBarrel = StringToEnumValue<EcalRecHit::Flags>(flagnamesbarrel);
473  auto const& flagnamesendcaps = cfg.getParameter<std::vector<std::string>>("recHitFlagsToBeExcludedEndcaps");
474  recHitsCfg.recHitFlagsToBeExcludedEndcaps = StringToEnumValue<EcalRecHit::Flags>(flagnamesendcaps);
475  auto const& severitynamesbarrel = cfg.getParameter<std::vector<std::string>>("recHitSeverityToBeExcludedBarrel");
476  recHitsCfg.recHitSeverityToBeExcludedBarrel =
477  StringToEnumValue<EcalSeverityLevel::SeverityLevel>(severitynamesbarrel);
478  auto const& severitynamesendcaps = cfg.getParameter<std::vector<std::string>>("recHitSeverityToBeExcludedEndcaps");
479  recHitsCfg.recHitSeverityToBeExcludedEndcaps =
480  StringToEnumValue<EcalSeverityLevel::SeverityLevel>(severitynamesendcaps);
481 
482  // isolation
484  .intRadiusHcal = cfg.getParameter<double>("intRadiusHcal"),
485  .etMinHcal = cfg.getParameter<double>("etMinHcal"),
486  .intRadiusEcalBarrel = cfg.getParameter<double>("intRadiusEcalBarrel"),
487  .intRadiusEcalEndcaps = cfg.getParameter<double>("intRadiusEcalEndcaps"),
488  .jurassicWidth = cfg.getParameter<double>("jurassicWidth"),
489  .etMinBarrel = cfg.getParameter<double>("etMinBarrel"),
490  .eMinBarrel = cfg.getParameter<double>("eMinBarrel"),
491  .etMinEndcaps = cfg.getParameter<double>("etMinEndcaps"),
492  .eMinEndcaps = cfg.getParameter<double>("eMinEndcaps"),
493  .vetoClustered = cfg.getParameter<bool>("vetoClustered"),
494  .useNumCrystals = cfg.getParameter<bool>("useNumCrystals")};
495 
496  // isolation
498  .ecaldrMax = pfECALClusIsolCfg.getParameter<double>("drMax"),
499  .ecaldrVetoBarrel = pfECALClusIsolCfg.getParameter<double>("drVetoBarrel"),
500  .ecaldrVetoEndcap = pfECALClusIsolCfg.getParameter<double>("drVetoEndcap"),
501  .ecaletaStripBarrel = pfECALClusIsolCfg.getParameter<double>("etaStripBarrel"),
502  .ecaletaStripEndcap = pfECALClusIsolCfg.getParameter<double>("etaStripEndcap"),
503  .ecalenergyBarrel = pfECALClusIsolCfg.getParameter<double>("energyBarrel"),
504  .ecalenergyEndcap = pfECALClusIsolCfg.getParameter<double>("energyEndcap"),
505  .useHF = pfHCALClusIsolCfg.getParameter<bool>("useHF"),
506  .hcaldrMax = pfHCALClusIsolCfg.getParameter<double>("drMax"),
507  .hcaldrVetoBarrel = pfHCALClusIsolCfg.getParameter<double>("drVetoBarrel"),
508  .hcaldrVetoEndcap = pfHCALClusIsolCfg.getParameter<double>("drVetoEndcap"),
509  .hcaletaStripBarrel = pfHCALClusIsolCfg.getParameter<double>("etaStripBarrel"),
510  .hcaletaStripEndcap = pfHCALClusIsolCfg.getParameter<double>("etaStripEndcap"),
511  .hcalenergyBarrel = pfHCALClusIsolCfg.getParameter<double>("energyBarrel"),
512  .hcalenergyEndcap = pfHCALClusIsolCfg.getParameter<double>("energyEndcap"),
513  .hcaluseEt = pfHCALClusIsolCfg.getParameter<bool>("useEt")};
514 
515  const RegressionHelper::Configuration regressionCfg{
516  .ecalRegressionWeightLabels = cfg.getParameter<std::vector<std::string>>("ecalRefinedRegressionWeightLabels"),
517  .ecalWeightsFromDB = cfg.getParameter<bool>("ecalWeightsFromDB"),
518  .ecalRegressionWeightFiles = cfg.getParameter<std::vector<std::string>>("ecalRefinedRegressionWeightFiles"),
519  .combinationRegressionWeightLabels =
520  cfg.getParameter<std::vector<std::string>>("combinationRegressionWeightLabels"),
521  .combinationWeightsFromDB = cfg.getParameter<bool>("combinationWeightsFromDB"),
522  .combinationRegressionWeightFiles =
523  cfg.getParameter<std::vector<std::string>>("combinationRegressionWeightFile")};
524 
525  // create algo
526  algo_ = std::make_unique<GsfElectronAlgo>(
527  inputCfg_,
528  strategyCfg_,
529  cutsCfg_,
530  hcalCfg_,
531  hcalCfgBc_,
532  isoCfg,
533  pfisoCfg,
534  recHitsCfg,
536  cfg.getParameter<std::string>("crackCorrectionFunction"), cfg, consumesCollector()),
537  regressionCfg,
538  cfg.getParameter<edm::ParameterSet>("trkIsol03Cfg"),
539  cfg.getParameter<edm::ParameterSet>("trkIsol04Cfg"),
540  cfg.getParameter<edm::ParameterSet>("trkIsolHEEP03Cfg"),
541  cfg.getParameter<edm::ParameterSet>("trkIsolHEEP04Cfg"),
542  consumesCollector());
543 
544  if (dnnPFidEnabled_) {
545  tfSessions_ = gcache->iElectronDNNEstimator->getSessions();
546  }
547 }
548 
550  for (auto session : tfSessions_) {
552  }
553 }
554 
556  if (!pset.exists("SeedConfiguration")) {
557  return;
558  }
559  edm::ParameterSet seedConfiguration = pset.getParameter<edm::ParameterSet>("SeedConfiguration");
560 
561  if (seedConfiguration.getParameter<bool>("applyHOverECut")) {
562  if ((hcalCfg_.hOverEConeSize != 0) &&
563  (hcalCfg_.hOverEConeSize != seedConfiguration.getParameter<double>("hOverEConeSize"))) {
564  edm::LogWarning("GsfElectronAlgo|InconsistentParameters")
565  << "The H/E cone size (" << hcalCfg_.hOverEConeSize << ") is different from ecal seeding ("
566  << seedConfiguration.getParameter<double>("hOverEConeSize") << ").";
567  }
568  if (cutsCfg_.maxHOverEBarrelCone < seedConfiguration.getParameter<double>("maxHOverEBarrel")) {
569  edm::LogWarning("GsfElectronAlgo|InconsistentParameters")
570  << "The max barrel cone H/E is lower than during ecal seeding.";
571  }
572  if (cutsCfg_.maxHOverEEndcapsCone < seedConfiguration.getParameter<double>("maxHOverEEndcaps")) {
573  edm::LogWarning("GsfElectronAlgo|InconsistentParameters")
574  << "The max endcaps cone H/E is lower than during ecal seeding.";
575  }
576  }
577 
578  if (cutsCfg_.minSCEtBarrel < seedConfiguration.getParameter<double>("SCEtCut")) {
579  edm::LogWarning("GsfElectronAlgo|InconsistentParameters")
580  << "The minimum super-cluster Et in barrel is lower than during ecal seeding.";
581  }
582  if (cutsCfg_.minSCEtEndcaps < seedConfiguration.getParameter<double>("SCEtCut")) {
583  edm::LogWarning("GsfElectronAlgo|InconsistentParameters")
584  << "The minimum super-cluster Et in endcaps is lower than during ecal seeding.";
585  }
586 }
587 
588 //=======================================================================================
589 // Ambiguity solving
590 //=======================================================================================
591 
593  edm::Event const& event,
594  bool ignoreNotPreselected) const {
595  // Getting required event data
596  auto const& beamspot = event.get(inputCfg_.beamSpotTag);
597  auto gsfPfRecTracks =
599  auto const& barrelRecHits = event.get(inputCfg_.barrelRecHitCollection);
600  auto const& endcapRecHits = event.get(inputCfg_.endcapRecHitCollection);
601 
602  if (strategyCfg_.ambSortingStrategy == 0) {
603  std::sort(electrons.begin(), electrons.end(), egamma::isBetterElectron);
604  } else if (strategyCfg_.ambSortingStrategy == 1) {
605  std::sort(electrons.begin(), electrons.end(), egamma::isInnermostElectron);
606  } else {
607  throw cms::Exception("GsfElectronAlgo|UnknownAmbiguitySortingStrategy")
608  << "value of strategyCfg_.ambSortingStrategy is : " << strategyCfg_.ambSortingStrategy;
609  }
610 
611  // init
612  for (auto& electron : electrons) {
613  electron.clearAmbiguousGsfTracks();
614  electron.setAmbiguous(false);
615  }
616 
617  // get ambiguous from GsfPfRecTracks
618  if (useGsfPfRecTracks_) {
619  for (auto& e1 : electrons) {
620  bool found = false;
621  for (auto const& gsfPfRecTrack : *gsfPfRecTracks) {
622  if (gsfPfRecTrack.gsfTrackRef() == e1.gsfTrack()) {
623  if (found) {
624  edm::LogWarning("GsfElectronAlgo") << "associated gsfPfRecTrack already found";
625  } else {
626  found = true;
627  for (auto const& duplicate : gsfPfRecTrack.convBremGsfPFRecTrackRef()) {
628  e1.addAmbiguousGsfTrack(duplicate->gsfTrackRef());
629  }
630  }
631  }
632  }
633  }
634  }
635  // or search overlapping clusters
636  else {
637  for (auto e1 = electrons.begin(); e1 != electrons.end(); ++e1) {
638  if (e1->ambiguous())
639  continue;
640  if (ignoreNotPreselected && !isPreselected(*e1))
641  continue;
642 
643  SuperClusterRef scRef1 = e1->superCluster();
644  CaloClusterPtr eleClu1 = e1->electronCluster();
645  LogDebug("GsfElectronAlgo") << "Blessing electron with E/P " << e1->eSuperClusterOverP() << ", cluster "
646  << scRef1.get() << " & track " << e1->gsfTrack().get();
647 
648  for (auto e2 = e1 + 1; e2 != electrons.end(); ++e2) {
649  if (e2->ambiguous())
650  continue;
651  if (ignoreNotPreselected && !isPreselected(*e2))
652  continue;
653 
654  SuperClusterRef scRef2 = e2->superCluster();
655  CaloClusterPtr eleClu2 = e2->electronCluster();
656 
657  // search if same cluster
658  bool sameCluster = false;
659  if (strategyCfg_.ambClustersOverlapStrategy == 0) {
660  sameCluster = (scRef1 == scRef2);
661  } else if (strategyCfg_.ambClustersOverlapStrategy == 1) {
662  float eMin = 1.;
663  float threshold = eMin * cosh(EleRelPoint(scRef1->position(), beamspot.position()).eta());
664  using egamma::sharedEnergy;
665  sameCluster = ((sharedEnergy(*eleClu1, *eleClu2, barrelRecHits, endcapRecHits) >= threshold) ||
666  (sharedEnergy(*scRef1->seed(), *eleClu2, barrelRecHits, endcapRecHits) >= threshold) ||
667  (sharedEnergy(*eleClu1, *scRef2->seed(), barrelRecHits, endcapRecHits) >= threshold) ||
668  (sharedEnergy(*scRef1->seed(), *scRef2->seed(), barrelRecHits, endcapRecHits) >= threshold));
669  } else {
670  throw cms::Exception("GsfElectronAlgo|UnknownAmbiguityClustersOverlapStrategy")
671  << "value of strategyCfg_.ambClustersOverlapStrategy is : " << strategyCfg_.ambClustersOverlapStrategy;
672  }
673 
674  // main instructions
675  if (sameCluster) {
676  LogDebug("GsfElectronAlgo") << "Discarding electron with E/P " << e2->eSuperClusterOverP() << ", cluster "
677  << scRef2.get() << " and track " << e2->gsfTrack().get();
678  e1->addAmbiguousGsfTrack(e2->gsfTrack());
679  e2->setAmbiguous(true);
680  } else if (e1->gsfTrack() == e2->gsfTrack()) {
681  edm::LogWarning("GsfElectronAlgo") << "Forgetting electron with E/P " << e2->eSuperClusterOverP()
682  << ", cluster " << scRef2.get() << " and track " << e2->gsfTrack().get();
683  e2->setAmbiguous(true);
684  }
685  }
686  }
687  }
688 }
689 
691  bool passCutBased = ele.passingCutBasedPreselection();
692  bool passPF = ele.passingPflowPreselection();
693  // it is worth nothing for gedGsfElectrons, this does nothing as its not set
694  // till GedGsfElectron finaliser, this is always false
695  bool passmva = ele.passingMvaPreselection();
696  if (!ele.ecalDrivenSeed()) {
697  if (ele.pt() > strategyCfg_.MaxElePtForOnlyMVA)
698  return passmva && passCutBased;
699  else
700  return passmva;
701  } else {
702  return passCutBased || passPF || passmva;
703  }
704 }
705 
706 // ------------ method called to produce the data ------------
708  // check configuration
711  auto seeds = event.getHandle(inputCfg_.seedsTag);
712  if (!seeds.isValid()) {
713  edm::LogWarning("GsfElectronAlgo|UnreachableSeedsProvenance")
714  << "Cannot check consistency of parameters with ecal seeding ones,"
715  << " because the original collection of seeds is not any more available.";
716  } else {
717  checkEcalSeedingParameters(edm::parameterSet(seeds.provenance()->stable(), event.processHistory()));
718  }
719  }
720 
721  auto electrons = algo_->completeElectrons(event, setup, globalCache());
723  const auto gsfMVAInputMap = matchWithPFCandidates(event.get(egmPFCandidateCollection_));
724  for (auto& el : electrons) {
725  el.setMvaInput(gsfMVAInputMap.find(el.gsfTrack())->second); // set Run2 MVA inputs
726  }
727  setMVAOutputs(electrons, globalCache(), event.get(inputCfg_.vtxCollectionTag), dnnPFidEnabled_, tfSessions_);
728  }
729 
730  // all electrons
731  logElectrons(electrons, event, "GsfElectronAlgo Info (before preselection)");
732  // preselection
733  if (strategyCfg_.applyPreselection) {
734  electrons.erase(
735  std::remove_if(electrons.begin(), electrons.end(), [this](auto const& ele) { return !isPreselected(ele); }),
736  electrons.end());
737  logElectrons(electrons, event, "GsfElectronAlgo Info (after preselection)");
738  }
739  // ambiguity
740  setAmbiguityData(electrons, event, strategyCfg_.ignoreNotPreselected);
741  if (strategyCfg_.applyAmbResolution) {
742  electrons.erase(std::remove_if(electrons.begin(), electrons.end(), std::mem_fn(&reco::GsfElectron::ambiguous)),
743  electrons.end());
744  logElectrons(electrons, event, "GsfElectronAlgo Info (after amb. solving)");
745  }
746  // go back to run2-like 2 effective depths if desired - depth 1 is the normal depth 1, depth 2 is the sum over the rest
747  if (hcalRun2EffDepth_) {
748  for (auto& ele : electrons)
749  ele.hcalToRun2EffDepth();
750  }
751  // final filling
752  event.emplace(electronPutToken_, std::move(electrons));
753 }
754 
void checkEcalSeedingParameters(edm::ParameterSet const &)
tuple cfg
Definition: looper.py:296
const edm::EDGetTokenT< reco::GsfPFRecTrackCollection > gsfPfRecTracksTag_
GsfElectronAlgo::Tokens inputCfg_
double pt() const final
transverse momentum
edm::EDGetTokenT< reco::PFCandidateCollection > egmPFCandidateCollection_
std::vector< std::string > ecalRegressionWeightLabels
void setAmbiguityData(reco::GsfElectronCollection &electrons, edm::Event const &event, bool ignoreNotPreselected=true) const
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
static void fillDescriptions(edm::ConfigurationDescriptions &)
const bool resetMvaValuesUsingPFCandidates_
bool exists(std::string const &parameterName) const
checks if a parameter exists
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
bool passingCutBasedPreselection() const
Definition: GsfElectron.h:763
bool passingMvaPreselection() const
Definition: GsfElectron.h:778
std::vector< GsfElectron > GsfElectronCollection
collection of GsfElectron objects
bool ambiguous() const
Definition: GsfElectron.h:765
#define LogTrace(id)
static std::string const input
Definition: EdmProvDump.cc:47
std::unique_ptr< const ElectronMVAEstimator > iElectronMVAEstimator
U second(std::pair< T, U > const &p)
bool isBetterElectron(reco::GsfElectron const &, reco::GsfElectron const &)
bool isInnermostElectron(reco::GsfElectron const &, reco::GsfElectron const &)
GsfElectronProducer(const edm::ParameterSet &, const GsfElectronAlgo::HeavyObjectCache *)
static std::unique_ptr< GsfElectronAlgo::HeavyObjectCache > initializeGlobalCache(const edm::ParameterSet &conf)
const edm::EDPutTokenT< reco::GsfElectronCollection > electronPutToken_
ParameterSet const & parameterSet(StableProvenance const &provenance, ProcessHistory const &history)
Definition: Provenance.cc:11
def move
Definition: eostools.py:511
bool closeSession(Session *&session)
Definition: TensorFlow.cc:198
bool get(ProductID const &oid, Handle< PROD > &result) const
Definition: Event.h:346
std::vector< tensorflow::Session * > tfSessions_
static edm::ParameterSetDescription pSetDescript()
T const * get() const
Returns C++ pointer to the item.
Definition: Ref.h:232
std::unique_ptr< const ElectronDNNEstimator > iElectronDNNEstimator
ParameterDescriptionBase * add(U const &iLabel, T const &value)
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
const GsfElectronAlgo::CutsConfiguration cutsCfg_
void add(std::string const &label, ParameterSetDescription const &psetDescription)
float sharedEnergy(reco::CaloCluster const &clu1, reco::CaloCluster const &clu2, EcalRecHitCollection const &barrelRecHits, EcalRecHitCollection const &endcapRecHits)
bool isPreselected(reco::GsfElectron const &ele) const
std::unique_ptr< GsfElectronAlgo > algo_
ElectronHcalHelper::Configuration hcalCfgBc_
std::unique_ptr< const SoftElectronMVAEstimator > sElectronMVAEstimator
#define get
Log< level::Warning, false > LogWarning
std::array< double, 4 > arrayHB
static void globalEndJob(GsfElectronAlgo::HeavyObjectCache const *)
void produce(edm::Event &event, const edm::EventSetup &setup) override
GsfElectronAlgo::StrategyConfiguration strategyCfg_
#define LogDebug(id)
bool ecalDrivenSeed() const
Definition: GsfElectron.h:158
bool passingPflowPreselection() const
Definition: GsfElectron.h:764
std::array< double, 7 > arrayHE