File indexing completed on 2021-04-12 23:30:34
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019 #include <memory>
0020
0021 #include "CommonTools/Utils/interface/KinematicTables.h"
0022 #include "CommonTools/Utils/interface/LazyConstructed.h"
0023 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0024 #include "DataFormats/EgammaReco/interface/ElectronSeed.h"
0025 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0026 #include "DataFormats/Math/interface/deltaR.h"
0027 #include "DataFormats/ParticleFlowReco/interface/GsfPFRecTrack.h"
0028 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0029 #include "DataFormats/ParticleFlowReco/interface/PFConversion.h"
0030 #include "DataFormats/ParticleFlowReco/interface/PFDisplacedTrackerVertex.h"
0031 #include "DataFormats/ParticleFlowReco/interface/PFDisplacedVertex.h"
0032 #include "DataFormats/ParticleFlowReco/interface/PFRecTrack.h"
0033 #include "DataFormats/ParticleFlowReco/interface/PFV0.h"
0034 #include "DataFormats/TrajectorySeed/interface/TrajectorySeed.h"
0035 #include "DataFormats/VertexReco/interface/Vertex.h"
0036 #include "FWCore/Framework/interface/ESHandle.h"
0037 #include "FWCore/Framework/interface/Event.h"
0038 #include "FWCore/Framework/interface/stream/EDProducer.h"
0039 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0040 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0041 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0042 #include "MagneticField/Engine/interface/MagneticField.h"
0043 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0044 #include "RecoParticleFlow/PFClusterTools/interface/ClusterClusterMapping.h"
0045 #include "RecoParticleFlow/PFClusterTools/interface/LinkByRecHit.h"
0046 #include "RecoParticleFlow/PFTracking/interface/ConvBremHeavyObjectCache.h"
0047 #include "RecoParticleFlow/PFTracking/interface/ConvBremPFTrackFinder.h"
0048 #include "RecoParticleFlow/PFTracking/interface/PFTrackAlgoTools.h"
0049 #include "RecoParticleFlow/PFTracking/interface/PFTrackTransformer.h"
0050 #include "TrackingTools/GsfTools/interface/MultiTrajectoryStateMode.h"
0051 #include "TrackingTools/GsfTools/interface/MultiTrajectoryStateTransform.h"
0052 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0053 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0054 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0055
0056 namespace {
0057
0058 constexpr float square(float x) { return x * x; };
0059
0060 }
0061
0062 class PFElecTkProducer final : public edm::stream::EDProducer<edm::GlobalCache<convbremhelpers::HeavyObjectCache> > {
0063 public:
0064
0065 explicit PFElecTkProducer(const edm::ParameterSet&, const convbremhelpers::HeavyObjectCache*);
0066
0067 static std::unique_ptr<convbremhelpers::HeavyObjectCache> initializeGlobalCache(const edm::ParameterSet& conf) {
0068 return std::make_unique<convbremhelpers::HeavyObjectCache>(conf);
0069 }
0070
0071 static void globalEndJob(convbremhelpers::HeavyObjectCache const*) {}
0072
0073 private:
0074 void beginRun(const edm::Run&, const edm::EventSetup&) override;
0075 void endRun(const edm::Run&, const edm::EventSetup&) override;
0076
0077
0078 void produce(edm::Event&, const edm::EventSetup&) override;
0079
0080 int findPfRef(const reco::PFRecTrackCollection& pfRTkColl,
0081 const reco::GsfTrack&,
0082 edm::soa::EtaPhiTableView trackEtaPhiTable);
0083
0084 bool applySelection(const reco::GsfTrack&);
0085
0086 bool resolveGsfTracks(const std::vector<reco::GsfPFRecTrack>& GsfPFVec,
0087 unsigned int ngsf,
0088 std::vector<unsigned int>& secondaries,
0089 const reco::PFClusterCollection& theEClus);
0090
0091 float minTangDist(const reco::GsfPFRecTrack& primGsf, const reco::GsfPFRecTrack& secGsf);
0092
0093 bool isSameEgSC(const reco::ElectronSeed& nSeed,
0094 const reco::ElectronSeed& iSeed,
0095 bool& bothGsfEcalDriven,
0096 float& SCEnergy);
0097
0098 bool isSharingEcalEnergyWithEgSC(const reco::GsfPFRecTrack& nGsfPFRecTrack,
0099 const reco::GsfPFRecTrack& iGsfPFRecTrack,
0100 const reco::ElectronSeed& nSeed,
0101 const reco::ElectronSeed& iSeed,
0102 const reco::PFClusterCollection& theEClus,
0103 bool& bothGsfTrackerDriven,
0104 bool& nEcalDriven,
0105 bool& iEcalDriven,
0106 float& nEnergy,
0107 float& iEnergy);
0108
0109 bool isInnerMost(const reco::GsfTrackRef& nGsfTrack, const reco::GsfTrackRef& iGsfTrack, bool& sameLayer);
0110
0111 bool isInnerMostWithLostHits(const reco::GsfTrackRef& nGsfTrack, const reco::GsfTrackRef& iGsfTrack, bool& sameLayer);
0112
0113 void createGsfPFRecTrackRef(const edm::OrphanHandle<reco::GsfPFRecTrackCollection>& gsfPfHandle,
0114 std::vector<reco::GsfPFRecTrack>& gsfPFRecTrackPrimary,
0115 const std::map<unsigned int, std::vector<reco::GsfPFRecTrack> >& MapPrimSec);
0116
0117
0118 reco::GsfPFRecTrack pftrack_;
0119 reco::GsfPFRecTrack secpftrack_;
0120 edm::ParameterSet conf_;
0121 edm::EDGetTokenT<reco::GsfTrackCollection> gsfTrackLabel_;
0122 edm::EDGetTokenT<reco::PFRecTrackCollection> pfTrackLabel_;
0123 edm::EDGetTokenT<reco::VertexCollection> primVtxLabel_;
0124 edm::EDGetTokenT<reco::PFClusterCollection> pfEcalClusters_;
0125 edm::EDGetTokenT<reco::PFDisplacedTrackerVertexCollection> pfNuclear_;
0126 edm::EDGetTokenT<reco::PFConversionCollection> pfConv_;
0127 edm::EDGetTokenT<reco::PFV0Collection> pfV0_;
0128 bool useNuclear_;
0129 bool useConversions_;
0130 bool useV0_;
0131 bool applyAngularGsfClean_;
0132 double detaCutGsfClean_;
0133 double dphiCutGsfClean_;
0134
0135 const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magFieldToken_;
0136 const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> tkerGeomToken_;
0137 const edm::ESGetToken<TransientTrackBuilder, TransientTrackRecord> transientTrackToken_;
0138
0139
0140 std::unique_ptr<PFTrackTransformer> pfTransformer_;
0141 MultiTrajectoryStateTransform mtsTransform_;
0142 std::unique_ptr<ConvBremPFTrackFinder> convBremFinder_;
0143
0144
0145 bool trajinev_;
0146 bool modemomentum_;
0147 bool applySel_;
0148 bool applyGsfClean_;
0149 bool useFifthStepForEcalDriven_;
0150 bool useFifthStepForTrackDriven_;
0151
0152 bool debugGsfCleaning_;
0153 double SCEne_;
0154 double detaGsfSC_;
0155 double dphiGsfSC_;
0156 double maxPtConvReco_;
0157
0158
0159 bool useConvBremFinder_;
0160
0161 double mvaConvBremFinderIDBarrelLowPt_;
0162 double mvaConvBremFinderIDBarrelHighPt_;
0163 double mvaConvBremFinderIDEndcapsLowPt_;
0164 double mvaConvBremFinderIDEndcapsHighPt_;
0165 std::string path_mvaWeightFileConvBremBarrelLowPt_;
0166 std::string path_mvaWeightFileConvBremBarrelHighPt_;
0167 std::string path_mvaWeightFileConvBremEndcapsLowPt_;
0168 std::string path_mvaWeightFileConvBremEndcapsHighPt_;
0169
0170
0171 std::vector<double> gsfInnerMomentumCache_;
0172 };
0173
0174 using namespace std;
0175 using namespace edm;
0176 using namespace reco;
0177
0178 PFElecTkProducer::PFElecTkProducer(const ParameterSet& iConfig, const convbremhelpers::HeavyObjectCache*)
0179 : conf_(iConfig),
0180 magFieldToken_(esConsumes<edm::Transition::BeginRun>()),
0181 tkerGeomToken_(esConsumes<edm::Transition::BeginRun>()),
0182 transientTrackToken_(esConsumes<edm::Transition::BeginRun>(edm::ESInputTag("", "TransientTrackBuilder"))) {
0183 gsfTrackLabel_ = consumes<reco::GsfTrackCollection>(iConfig.getParameter<InputTag>("GsfTrackModuleLabel"));
0184
0185 pfTrackLabel_ = consumes<reco::PFRecTrackCollection>(iConfig.getParameter<InputTag>("PFRecTrackLabel"));
0186
0187 primVtxLabel_ = consumes<reco::VertexCollection>(iConfig.getParameter<InputTag>("PrimaryVertexLabel"));
0188
0189 pfEcalClusters_ = consumes<reco::PFClusterCollection>(iConfig.getParameter<InputTag>("PFEcalClusters"));
0190
0191 pfNuclear_ = consumes<reco::PFDisplacedTrackerVertexCollection>(iConfig.getParameter<InputTag>("PFNuclear"));
0192
0193 pfConv_ = consumes<reco::PFConversionCollection>(iConfig.getParameter<InputTag>("PFConversions"));
0194
0195 pfV0_ = consumes<reco::PFV0Collection>(iConfig.getParameter<InputTag>("PFV0"));
0196
0197 useNuclear_ = iConfig.getParameter<bool>("useNuclear");
0198 useConversions_ = iConfig.getParameter<bool>("useConversions");
0199 useV0_ = iConfig.getParameter<bool>("useV0");
0200 debugGsfCleaning_ = iConfig.getParameter<bool>("debugGsfCleaning");
0201
0202 produces<GsfPFRecTrackCollection>();
0203 produces<GsfPFRecTrackCollection>("Secondary").setBranchAlias("secondary");
0204
0205 trajinev_ = iConfig.getParameter<bool>("TrajInEvents");
0206 modemomentum_ = iConfig.getParameter<bool>("ModeMomentum");
0207 applySel_ = iConfig.getParameter<bool>("applyEGSelection");
0208 applyGsfClean_ = iConfig.getParameter<bool>("applyGsfTrackCleaning");
0209 applyAngularGsfClean_ = iConfig.getParameter<bool>("applyAlsoGsfAngularCleaning");
0210 detaCutGsfClean_ = iConfig.getParameter<double>("maxDEtaGsfAngularCleaning");
0211 dphiCutGsfClean_ = iConfig.getParameter<double>("maxDPhiBremTangGsfAngularCleaning");
0212 useFifthStepForTrackDriven_ = iConfig.getParameter<bool>("useFifthStepForTrackerDrivenGsf");
0213 useFifthStepForEcalDriven_ = iConfig.getParameter<bool>("useFifthStepForEcalDrivenGsf");
0214 maxPtConvReco_ = iConfig.getParameter<double>("MaxConvBremRecoPT");
0215 detaGsfSC_ = iConfig.getParameter<double>("MinDEtaGsfSC");
0216 dphiGsfSC_ = iConfig.getParameter<double>("MinDPhiGsfSC");
0217 SCEne_ = iConfig.getParameter<double>("MinSCEnergy");
0218
0219
0220 useConvBremFinder_ = iConfig.getParameter<bool>("useConvBremFinder");
0221
0222 mvaConvBremFinderIDBarrelLowPt_ = iConfig.getParameter<double>("pf_convBremFinderID_mvaCutBarrelLowPt");
0223 mvaConvBremFinderIDBarrelHighPt_ = iConfig.getParameter<double>("pf_convBremFinderID_mvaCutBarrelHighPt");
0224 mvaConvBremFinderIDEndcapsLowPt_ = iConfig.getParameter<double>("pf_convBremFinderID_mvaCutEndcapsLowPt");
0225 mvaConvBremFinderIDEndcapsHighPt_ = iConfig.getParameter<double>("pf_convBremFinderID_mvaCutEndcapsHighPt");
0226 }
0227
0228
0229
0230
0231
0232
0233 void PFElecTkProducer::produce(Event& iEvent, const EventSetup& iSetup) {
0234
0235 auto gsfPFRecTrackCollection = std::make_unique<GsfPFRecTrackCollection>();
0236
0237 auto gsfPFRecTrackCollectionSecondary = std::make_unique<GsfPFRecTrackCollection>();
0238
0239
0240 Handle<GsfTrackCollection> gsftrackscoll;
0241 iEvent.getByToken(gsfTrackLabel_, gsftrackscoll);
0242
0243
0244 Handle<vector<Trajectory> > TrajectoryCollection;
0245
0246
0247 Handle<PFRecTrackCollection> thePfRecTrackCollection;
0248 iEvent.getByToken(pfTrackLabel_, thePfRecTrackCollection);
0249
0250
0251
0252 auto trackEtaPhiTable = makeLazy<edm::soa::EtaPhiTable>(*thePfRecTrackCollection);
0253
0254
0255 Handle<PFClusterCollection> theECPfClustCollection;
0256 iEvent.getByToken(pfEcalClusters_, theECPfClustCollection);
0257 const PFClusterCollection& theEcalClusters = *(theECPfClustCollection.product());
0258
0259
0260 Handle<reco::VertexCollection> thePrimaryVertexColl;
0261 iEvent.getByToken(primVtxLabel_, thePrimaryVertexColl);
0262
0263
0264 Handle<reco::PFDisplacedTrackerVertexCollection> pfNuclears;
0265 if (useNuclear_)
0266 iEvent.getByToken(pfNuclear_, pfNuclears);
0267
0268
0269 Handle<reco::PFConversionCollection> pfConversions;
0270 if (useConversions_)
0271 iEvent.getByToken(pfConv_, pfConversions);
0272
0273
0274 Handle<reco::PFV0Collection> pfV0;
0275 if (useV0_)
0276 iEvent.getByToken(pfV0_, pfV0);
0277
0278 GsfTrackCollection gsftracks = *(gsftrackscoll.product());
0279 vector<Trajectory> tjvec(0);
0280 if (trajinev_) {
0281 iEvent.getByToken(gsfTrackLabel_, TrajectoryCollection);
0282
0283 tjvec = *(TrajectoryCollection.product());
0284 }
0285
0286 vector<reco::GsfPFRecTrack> selGsfPFRecTracks;
0287 vector<reco::GsfPFRecTrack> primaryGsfPFRecTracks;
0288 std::map<unsigned int, std::vector<reco::GsfPFRecTrack> > GsfPFMap;
0289
0290 for (unsigned igsf = 0; igsf < gsftracks.size(); igsf++) {
0291 GsfTrackRef trackRef(gsftrackscoll, igsf);
0292 gsfInnerMomentumCache_.push_back(trackRef->pMode());
0293 TrajectoryStateOnSurface i_inTSOS = mtsTransform_.innerStateOnSurface((*trackRef));
0294 GlobalVector i_innMom;
0295 if (i_inTSOS.isValid()) {
0296 multiTrajectoryStateMode::momentumFromModeCartesian(i_inTSOS, i_innMom);
0297 gsfInnerMomentumCache_.back() = i_innMom.mag();
0298 }
0299 }
0300
0301 for (unsigned int igsf = 0; igsf < gsftracks.size(); igsf++) {
0302 GsfTrackRef trackRef(gsftrackscoll, igsf);
0303
0304 int kf_ind = findPfRef(*thePfRecTrackCollection, gsftracks[igsf], trackEtaPhiTable.value());
0305
0306 if (kf_ind >= 0) {
0307 PFRecTrackRef kf_ref(thePfRecTrackCollection, kf_ind);
0308
0309
0310 if (useFifthStepForEcalDriven_ == false || useFifthStepForTrackDriven_ == false) {
0311 bool isFifthStepTrack = PFTrackAlgoTools::isFifthStep(kf_ref->trackRef()->algo());
0312 bool isEcalDriven = true;
0313 bool isTrackerDriven = true;
0314
0315 if (trackRef->seedRef().get() == nullptr) {
0316 isEcalDriven = false;
0317 isTrackerDriven = false;
0318 } else {
0319 auto const& SeedFromRef = static_cast<ElectronSeed const&>(*(trackRef->extra()->seedRef()));
0320 if (SeedFromRef.caloCluster().isNull())
0321 isEcalDriven = false;
0322 if (SeedFromRef.ctfTrack().isNull())
0323 isTrackerDriven = false;
0324 }
0325
0326 if (isFifthStepTrack && isEcalDriven && isTrackerDriven == false && useFifthStepForEcalDriven_ == false) {
0327 continue;
0328 }
0329
0330 if (isFifthStepTrack && isTrackerDriven && isEcalDriven == false && useFifthStepForTrackDriven_ == false) {
0331 continue;
0332 }
0333
0334 if (isFifthStepTrack && isTrackerDriven && isEcalDriven && useFifthStepForTrackDriven_ == false &&
0335 useFifthStepForEcalDriven_ == false) {
0336 continue;
0337 }
0338 }
0339
0340 pftrack_ = GsfPFRecTrack(gsftracks[igsf].charge(), reco::PFRecTrack::GSF, igsf, trackRef, kf_ref);
0341 } else {
0342 PFRecTrackRef dummyRef;
0343 pftrack_ = GsfPFRecTrack(gsftracks[igsf].charge(), reco::PFRecTrack::GSF, igsf, trackRef, dummyRef);
0344 }
0345
0346 bool validgsfbrem = false;
0347 if (trajinev_) {
0348 validgsfbrem = pfTransformer_->addPointsAndBrems(pftrack_, gsftracks[igsf], tjvec[igsf], modemomentum_);
0349 } else {
0350 validgsfbrem = pfTransformer_->addPointsAndBrems(pftrack_, gsftracks[igsf], mtsTransform_);
0351 }
0352
0353 bool passSel = true;
0354 if (applySel_)
0355 passSel = applySelection(gsftracks[igsf]);
0356
0357 if (validgsfbrem && passSel)
0358 selGsfPFRecTracks.push_back(pftrack_);
0359 }
0360
0361 unsigned int count_primary = 0;
0362 if (!selGsfPFRecTracks.empty()) {
0363 for (unsigned int ipfgsf = 0; ipfgsf < selGsfPFRecTracks.size(); ipfgsf++) {
0364 vector<unsigned int> secondaries(0);
0365 secondaries.clear();
0366 bool keepGsf = true;
0367
0368 if (applyGsfClean_) {
0369 keepGsf = resolveGsfTracks(selGsfPFRecTracks, ipfgsf, secondaries, theEcalClusters);
0370 }
0371
0372
0373 if (keepGsf == true) {
0374
0375 if (convBremFinder_->foundConvBremPFRecTrack(thePfRecTrackCollection,
0376 thePrimaryVertexColl,
0377 pfNuclears,
0378 pfConversions,
0379 pfV0,
0380 globalCache(),
0381 useNuclear_,
0382 useConversions_,
0383 useV0_,
0384 theEcalClusters,
0385 selGsfPFRecTracks[ipfgsf])) {
0386 const vector<PFRecTrackRef>& convBremPFRecTracks(convBremFinder_->getConvBremPFRecTracks());
0387 for (unsigned int ii = 0; ii < convBremPFRecTracks.size(); ii++) {
0388 selGsfPFRecTracks[ipfgsf].addConvBremPFRecTrackRef(convBremPFRecTracks[ii]);
0389 }
0390 }
0391
0392
0393
0394 primaryGsfPFRecTracks.push_back(selGsfPFRecTracks[ipfgsf]);
0395
0396
0397
0398
0399 unsigned int primGsfIndex = selGsfPFRecTracks[ipfgsf].trackId();
0400 vector<reco::GsfPFRecTrack> trueGsfPFRecTracks;
0401 if (!secondaries.empty()) {
0402
0403 for (unsigned int isecpfgsf = 0; isecpfgsf < secondaries.size(); isecpfgsf++) {
0404 PFRecTrackRef refsecKF = selGsfPFRecTracks[(secondaries[isecpfgsf])].kfPFRecTrackRef();
0405
0406 unsigned int secGsfIndex = selGsfPFRecTracks[(secondaries[isecpfgsf])].trackId();
0407 GsfTrackRef secGsfRef = selGsfPFRecTracks[(secondaries[isecpfgsf])].gsfTrackRef();
0408
0409 if (refsecKF.isNonnull()) {
0410
0411 secpftrack_ = GsfPFRecTrack(
0412 gsftracks[secGsfIndex].charge(), reco::PFRecTrack::GSF, primGsfIndex, secGsfRef, refsecKF);
0413 } else {
0414 PFRecTrackRef dummyRef;
0415
0416 secpftrack_ = GsfPFRecTrack(
0417 gsftracks[secGsfIndex].charge(), reco::PFRecTrack::GSF, primGsfIndex, secGsfRef, dummyRef);
0418 }
0419
0420 bool validgsfbrem = false;
0421 if (trajinev_) {
0422 validgsfbrem = pfTransformer_->addPointsAndBrems(
0423 secpftrack_, gsftracks[secGsfIndex], tjvec[secGsfIndex], modemomentum_);
0424 } else {
0425 validgsfbrem = pfTransformer_->addPointsAndBrems(secpftrack_, gsftracks[secGsfIndex], mtsTransform_);
0426 }
0427
0428 if (validgsfbrem) {
0429 gsfPFRecTrackCollectionSecondary->push_back(secpftrack_);
0430 trueGsfPFRecTracks.push_back(secpftrack_);
0431 }
0432 }
0433 }
0434 GsfPFMap.insert(pair<unsigned int, std::vector<reco::GsfPFRecTrack> >(count_primary, trueGsfPFRecTracks));
0435 trueGsfPFRecTracks.clear();
0436 count_primary++;
0437 }
0438 }
0439 }
0440
0441 const edm::OrphanHandle<GsfPFRecTrackCollection> gsfPfRefProd =
0442 iEvent.put(std::move(gsfPFRecTrackCollectionSecondary), "Secondary");
0443
0444
0445 createGsfPFRecTrackRef(gsfPfRefProd, primaryGsfPFRecTracks, GsfPFMap);
0446
0447 for (unsigned int iGSF = 0; iGSF < primaryGsfPFRecTracks.size(); iGSF++) {
0448 gsfPFRecTrackCollection->push_back(primaryGsfPFRecTracks[iGSF]);
0449 }
0450 iEvent.put(std::move(gsfPFRecTrackCollection));
0451
0452 selGsfPFRecTracks.clear();
0453 GsfPFMap.clear();
0454 primaryGsfPFRecTracks.clear();
0455
0456 std::vector<double>().swap(gsfInnerMomentumCache_);
0457 }
0458
0459
0460 void PFElecTkProducer::createGsfPFRecTrackRef(
0461 const edm::OrphanHandle<reco::GsfPFRecTrackCollection>& gsfPfHandle,
0462 std::vector<reco::GsfPFRecTrack>& gsfPFRecTrackPrimary,
0463 const std::map<unsigned int, std::vector<reco::GsfPFRecTrack> >& MapPrimSec) {
0464 unsigned int cgsf = 0;
0465 unsigned int csecgsf = 0;
0466 for (std::map<unsigned int, std::vector<reco::GsfPFRecTrack> >::const_iterator igsf = MapPrimSec.begin();
0467 igsf != MapPrimSec.end();
0468 igsf++, cgsf++) {
0469 vector<reco::GsfPFRecTrack> SecGsfPF = igsf->second;
0470 for (unsigned int iSecGsf = 0; iSecGsf < SecGsfPF.size(); iSecGsf++) {
0471 edm::Ref<reco::GsfPFRecTrackCollection> refgprt(gsfPfHandle, csecgsf);
0472 gsfPFRecTrackPrimary[cgsf].addConvBremGsfPFRecTrackRef(refgprt);
0473 ++csecgsf;
0474 }
0475 }
0476
0477 return;
0478 }
0479
0480 int PFElecTkProducer::findPfRef(const reco::PFRecTrackCollection& pfRTkColl,
0481 const reco::GsfTrack& gsftk,
0482 edm::soa::EtaPhiTableView trackEtaPhiTable) {
0483 if (gsftk.seedRef().get() == nullptr) {
0484 return -1;
0485 }
0486
0487
0488 constexpr float maxDR2 = square(0.05f);
0489
0490
0491 float gsftkEta = gsftk.eta();
0492 float gsftkPhi = gsftk.phi();
0493 auto const& gsftkHits = gsftk.seedRef()->recHits();
0494
0495 auto const& electronSeedFromRef = static_cast<ElectronSeed const&>(*(gsftk.extra()->seedRef()));
0496
0497 if (electronSeedFromRef.ctfTrack().isNull()) {
0498 unsigned int i_pf = 0;
0499 int ibest = -1;
0500 unsigned int ish_max = 0;
0501 float dr2_min = square(1000.f);
0502
0503
0504
0505 for (auto const& pft : trackEtaPhiTable) {
0506 unsigned int ish = 0;
0507
0508 using namespace edm::soa::col;
0509 const float dr2 = reco::deltaR2(pft.get<Eta>(), pft.get<Phi>(), gsftkEta, gsftkPhi);
0510
0511 if (dr2 <= maxDR2) {
0512 for (auto const& hhit : pfRTkColl[i_pf].trackRef()->recHits()) {
0513 if (!hhit->isValid())
0514 continue;
0515 for (auto const& hit : gsftkHits) {
0516 if (hit.isValid() && hhit->sharesInput(&hit, TrackingRecHit::all))
0517 ish++;
0518 }
0519 }
0520
0521 if ((ish > ish_max) || ((ish == ish_max) && (dr2 < dr2_min))) {
0522 ish_max = ish;
0523 dr2_min = dr2;
0524 ibest = i_pf;
0525 }
0526 }
0527
0528 i_pf++;
0529 }
0530
0531 return ((ish_max == 0) || (dr2_min > maxDR2)) ? -1 : ibest;
0532 } else {
0533
0534
0535 unsigned int i_pf = 0;
0536
0537 for (auto const& pft : pfRTkColl) {
0538
0539 if (pft.trackRef() == electronSeedFromRef.ctfTrack()) {
0540 return i_pf;
0541 }
0542 i_pf++;
0543 }
0544 }
0545 return -1;
0546 }
0547
0548
0549 bool PFElecTkProducer::applySelection(const reco::GsfTrack& gsftk) {
0550 if (gsftk.seedRef().get() == nullptr)
0551 return false;
0552 auto const& ElSeedFromRef = static_cast<ElectronSeed const&>(*(gsftk.extra()->seedRef()));
0553
0554 bool passCut = false;
0555 if (ElSeedFromRef.ctfTrack().isNull()) {
0556 if (ElSeedFromRef.caloCluster().isNull())
0557 return passCut;
0558 auto const* scRef = static_cast<SuperCluster const*>(ElSeedFromRef.caloCluster().get());
0559
0560 if (scRef) {
0561 float caloEne = scRef->energy();
0562 float feta = fabs(scRef->eta() - gsftk.etaMode());
0563 float fphi = fabs(scRef->phi() - gsftk.phiMode());
0564 if (fphi > TMath::Pi())
0565 fphi -= TMath::TwoPi();
0566 if (caloEne > SCEne_ && feta < detaGsfSC_ && fabs(fphi) < dphiGsfSC_)
0567 passCut = true;
0568 }
0569 } else {
0570
0571 passCut = true;
0572 }
0573 return passCut;
0574 }
0575 bool PFElecTkProducer::resolveGsfTracks(const vector<reco::GsfPFRecTrack>& GsfPFVec,
0576 unsigned int ngsf,
0577 vector<unsigned int>& secondaries,
0578 const reco::PFClusterCollection& theEClus) {
0579 bool debugCleaning = debugGsfCleaning_;
0580 bool n_keepGsf = true;
0581
0582 const reco::GsfTrackRef& nGsfTrack = GsfPFVec[ngsf].gsfTrackRef();
0583
0584 if (nGsfTrack->seedRef().get() == nullptr)
0585 return false;
0586 auto const& nElSeedFromRef = static_cast<ElectronSeed const&>(*(nGsfTrack->extra()->seedRef()));
0587
0588
0589
0590
0591
0592
0593
0594
0595
0596
0597 float nPin = gsfInnerMomentumCache_[nGsfTrack.key()];
0598
0599 float neta = nGsfTrack->innerMomentum().eta();
0600 float nphi = nGsfTrack->innerMomentum().phi();
0601
0602 if (debugCleaning)
0603 cout << " PFElecTkProducer:: considering track " << nGsfTrack->pt() << " eta,phi " << nGsfTrack->eta() << ", "
0604 << nGsfTrack->phi() << endl;
0605
0606 for (unsigned int igsf = 0; igsf < GsfPFVec.size(); igsf++) {
0607 if (igsf != ngsf) {
0608 reco::GsfTrackRef iGsfTrack = GsfPFVec[igsf].gsfTrackRef();
0609
0610 if (debugCleaning)
0611 cout << " PFElecTkProducer:: and comparing with track " << iGsfTrack->pt() << " eta,phi " << iGsfTrack->eta()
0612 << ", " << iGsfTrack->phi() << endl;
0613
0614 float ieta = iGsfTrack->innerMomentum().eta();
0615 float iphi = iGsfTrack->innerMomentum().phi();
0616 float feta = fabs(neta - ieta);
0617 float fphi = fabs(nphi - iphi);
0618 if (fphi > TMath::Pi())
0619 fphi -= TMath::TwoPi();
0620
0621
0622 if (feta < 0.5 && fabs(fphi) < 1.0) {
0623 if (debugCleaning)
0624 cout << " Entering angular superloose preselection " << endl;
0625
0626
0627
0628
0629
0630
0631
0632
0633
0634
0635 float iPin = gsfInnerMomentumCache_[iGsfTrack.key()];
0636
0637 if (iGsfTrack->seedRef().get() == nullptr)
0638 continue;
0639 auto const& iElSeedFromRef = static_cast<ElectronSeed const&>(*(iGsfTrack->extra()->seedRef()));
0640
0641 float SCEnergy = -1.;
0642
0643 bool areBothGsfEcalDriven = false;
0644 ;
0645 bool isSameSC = isSameEgSC(nElSeedFromRef, iElSeedFromRef, areBothGsfEcalDriven, SCEnergy);
0646
0647
0648 if (areBothGsfEcalDriven) {
0649 if (isSameSC) {
0650 float nEP = SCEnergy / nPin;
0651 float iEP = SCEnergy / iPin;
0652 if (debugCleaning)
0653 cout << " Entering SAME supercluster case "
0654 << " nEP " << nEP << " iEP " << iEP << endl;
0655
0656
0657
0658
0659
0660 bool isSameLayer = false;
0661 bool iGsfInnermostWithLostHits = isInnerMostWithLostHits(nGsfTrack, iGsfTrack, isSameLayer);
0662
0663 if (debugCleaning)
0664 cout << " iGsf is InnerMostWithLostHits " << iGsfInnermostWithLostHits << " isSameLayer " << isSameLayer
0665 << endl;
0666
0667 if (iGsfInnermostWithLostHits) {
0668 n_keepGsf = false;
0669 return n_keepGsf;
0670 } else if (isSameLayer) {
0671 if (fabs(iEP - 1) < fabs(nEP - 1)) {
0672 n_keepGsf = false;
0673 return n_keepGsf;
0674 } else {
0675 secondaries.push_back(igsf);
0676 }
0677 } else {
0678
0679 secondaries.push_back(igsf);
0680 }
0681 }
0682 } else {
0683
0684 float minBremDphi = minTangDist(GsfPFVec[ngsf], GsfPFVec[igsf]);
0685 float nETot = 0.;
0686 float iETot = 0.;
0687 bool isBothGsfTrackerDriven = false;
0688 bool nEcalDriven = false;
0689 bool iEcalDriven = false;
0690 bool isSameScEgPf = isSharingEcalEnergyWithEgSC(GsfPFVec[ngsf],
0691 GsfPFVec[igsf],
0692 nElSeedFromRef,
0693 iElSeedFromRef,
0694 theEClus,
0695 isBothGsfTrackerDriven,
0696 nEcalDriven,
0697 iEcalDriven,
0698 nETot,
0699 iETot);
0700
0701
0702 bool isSameLayer = false;
0703 bool iGsfInnermostWithLostHits = isInnerMostWithLostHits(nGsfTrack, iGsfTrack, isSameLayer);
0704
0705 if (isSameScEgPf) {
0706
0707
0708 if (debugCleaning) {
0709 cout << " Sharing ECAL energy passed "
0710 << " nEtot " << nETot << " iEtot " << iETot << endl;
0711 if (isBothGsfTrackerDriven)
0712 cout << " Both Track are trackerDriven " << endl;
0713 }
0714
0715
0716 if (iGsfInnermostWithLostHits) {
0717 n_keepGsf = false;
0718 return n_keepGsf;
0719 } else if (isSameLayer) {
0720
0721
0722
0723
0724
0725
0726
0727 if (isBothGsfTrackerDriven == false) {
0728
0729 if (iEcalDriven) {
0730 n_keepGsf = false;
0731 return n_keepGsf;
0732 } else {
0733 secondaries.push_back(igsf);
0734 }
0735 } else {
0736
0737
0738
0739 float ETot = -1;
0740 if (nETot != iETot) {
0741 if (nETot > iETot)
0742 ETot = nETot;
0743 else
0744 ETot = iETot;
0745 } else {
0746 ETot = nETot;
0747 }
0748 float nEP = ETot / nPin;
0749 float iEP = ETot / iPin;
0750
0751 if (debugCleaning)
0752 cout << " nETot " << nETot << " iETot " << iETot << " ETot " << ETot << endl
0753 << " nPin " << nPin << " iPin " << iPin << " nEP " << nEP << " iEP " << iEP << endl;
0754
0755 if (fabs(iEP - 1) < fabs(nEP - 1)) {
0756 n_keepGsf = false;
0757 return n_keepGsf;
0758 } else {
0759 secondaries.push_back(igsf);
0760 }
0761 }
0762 } else {
0763 secondaries.push_back(igsf);
0764 }
0765 } else if (feta < detaCutGsfClean_ && minBremDphi < dphiCutGsfClean_) {
0766
0767 bool secPushedBack = false;
0768 if (nEcalDriven == false && nETot == 0.) {
0769 n_keepGsf = false;
0770 return n_keepGsf;
0771 } else if (iEcalDriven == false && iETot == 0.) {
0772 secondaries.push_back(igsf);
0773 secPushedBack = true;
0774 }
0775 if (debugCleaning)
0776 cout << " Close Tracks "
0777 << " feta " << feta << " fabs(fphi) " << fabs(fphi) << " minBremDphi " << minBremDphi << " nETot "
0778 << nETot << " iETot " << iETot << " nLostHits " << nGsfTrack->missingInnerHits() << " iLostHits "
0779 << iGsfTrack->missingInnerHits() << endl;
0780
0781
0782 if (applyAngularGsfClean_) {
0783 if (iGsfInnermostWithLostHits) {
0784 n_keepGsf = false;
0785 return n_keepGsf;
0786 } else if (isSameLayer == false) {
0787 if (secPushedBack == false)
0788 secondaries.push_back(igsf);
0789 }
0790 }
0791 } else if (feta < 0.1 && minBremDphi < 0.2) {
0792
0793
0794 if (debugCleaning)
0795 cout << " Close Tracks and failed all the conditions "
0796 << " feta " << feta << " fabs(fphi) " << fabs(fphi) << " minBremDphi " << minBremDphi << " nETot "
0797 << nETot << " iETot " << iETot << " nLostHits " << nGsfTrack->missingInnerHits() << " iLostHits "
0798 << iGsfTrack->missingInnerHits() << endl;
0799
0800 if (nEcalDriven == false && nETot == 0.) {
0801 n_keepGsf = false;
0802 return n_keepGsf;
0803 }
0804
0805 }
0806 }
0807 }
0808 }
0809 }
0810
0811 return n_keepGsf;
0812 }
0813 float PFElecTkProducer::minTangDist(const reco::GsfPFRecTrack& primGsf, const reco::GsfPFRecTrack& secGsf) {
0814 float minDphi = 1000.;
0815
0816 std::vector<reco::PFBrem> primPFBrem = primGsf.PFRecBrem();
0817 std::vector<reco::PFBrem> secPFBrem = secGsf.PFRecBrem();
0818
0819 unsigned int cbrem = 0;
0820 for (unsigned isbrem = 0; isbrem < secPFBrem.size(); isbrem++) {
0821 if (secPFBrem[isbrem].indTrajPoint() == 99)
0822 continue;
0823 const reco::PFTrajectoryPoint& atSecECAL =
0824 secPFBrem[isbrem].extrapolatedPoint(reco::PFTrajectoryPoint::ECALEntrance);
0825 if (!atSecECAL.isValid())
0826 continue;
0827 float secPhi = atSecECAL.positionREP().Phi();
0828
0829 unsigned int sbrem = 0;
0830 for (unsigned ipbrem = 0; ipbrem < primPFBrem.size(); ipbrem++) {
0831 if (primPFBrem[ipbrem].indTrajPoint() == 99)
0832 continue;
0833 const reco::PFTrajectoryPoint& atPrimECAL =
0834 primPFBrem[ipbrem].extrapolatedPoint(reco::PFTrajectoryPoint::ECALEntrance);
0835 if (!atPrimECAL.isValid())
0836 continue;
0837 sbrem++;
0838 if (sbrem <= 3) {
0839 float primPhi = atPrimECAL.positionREP().Phi();
0840
0841 float dphi = fabs(primPhi - secPhi);
0842 if (dphi > TMath::Pi())
0843 dphi -= TMath::TwoPi();
0844 if (fabs(dphi) < minDphi) {
0845 minDphi = fabs(dphi);
0846 }
0847 }
0848 }
0849
0850 cbrem++;
0851 if (cbrem == 3)
0852 break;
0853 }
0854 return minDphi;
0855 }
0856 bool PFElecTkProducer::isSameEgSC(const reco::ElectronSeed& nSeed,
0857 const reco::ElectronSeed& iSeed,
0858 bool& bothGsfEcalDriven,
0859 float& SCEnergy) {
0860 bool isSameSC = false;
0861
0862 if (nSeed.caloCluster().isNonnull() && iSeed.caloCluster().isNonnull()) {
0863 auto const* nscRef = static_cast<SuperCluster const*>(nSeed.caloCluster().get());
0864 auto const* iscRef = static_cast<SuperCluster const*>(iSeed.caloCluster().get());
0865
0866 if (nscRef && iscRef) {
0867 bothGsfEcalDriven = true;
0868 if (nscRef == iscRef) {
0869 isSameSC = true;
0870
0871 SCEnergy = nscRef->energy();
0872 }
0873 }
0874 }
0875 return isSameSC;
0876 }
0877 bool PFElecTkProducer::isSharingEcalEnergyWithEgSC(const reco::GsfPFRecTrack& nGsfPFRecTrack,
0878 const reco::GsfPFRecTrack& iGsfPFRecTrack,
0879 const reco::ElectronSeed& nSeed,
0880 const reco::ElectronSeed& iSeed,
0881 const reco::PFClusterCollection& theEClus,
0882 bool& bothGsfTrackerDriven,
0883 bool& nEcalDriven,
0884 bool& iEcalDriven,
0885 float& nEnergy,
0886 float& iEnergy) {
0887 bool isSharingEnergy = false;
0888
0889
0890 bool oneEcalDriven = true;
0891 SuperCluster const* scRef = nullptr;
0892 GsfPFRecTrack gsfPfTrack;
0893
0894 if (nSeed.caloCluster().isNonnull()) {
0895 scRef = static_cast<SuperCluster const*>(nSeed.caloCluster().get());
0896 assert(scRef);
0897 nEnergy = scRef->energy();
0898 nEcalDriven = true;
0899 gsfPfTrack = iGsfPFRecTrack;
0900 } else if (iSeed.caloCluster().isNonnull()) {
0901 scRef = static_cast<SuperCluster const*>(iSeed.caloCluster().get());
0902 assert(scRef);
0903 iEnergy = scRef->energy();
0904 iEcalDriven = true;
0905 gsfPfTrack = nGsfPFRecTrack;
0906 } else {
0907 oneEcalDriven = false;
0908 }
0909
0910 if (oneEcalDriven) {
0911
0912
0913 vector<PFCluster> vecPFClusters;
0914 vecPFClusters.clear();
0915
0916 for (PFClusterCollection::const_iterator clus = theEClus.begin(); clus != theEClus.end(); clus++) {
0917 PFCluster clust = *clus;
0918 clust.calculatePositionREP();
0919
0920 float deta = fabs(scRef->position().eta() - clust.position().eta());
0921 float dphi = fabs(scRef->position().phi() - clust.position().phi());
0922 if (dphi > TMath::Pi())
0923 dphi -= TMath::TwoPi();
0924
0925
0926
0927 if (deta < 0.5 && fabs(dphi) < 1.0) {
0928 double distGsf = gsfPfTrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALShowerMax).isValid()
0929 ? LinkByRecHit::testTrackAndClusterByRecHit(gsfPfTrack, clust)
0930 : -1.;
0931
0932 if (distGsf > 0.) {
0933 if (nEcalDriven)
0934 iEnergy += clust.energy();
0935 else
0936 nEnergy += clust.energy();
0937 vecPFClusters.push_back(clust);
0938 }
0939
0940 else {
0941 vector<PFBrem> primPFBrem = gsfPfTrack.PFRecBrem();
0942 for (unsigned ipbrem = 0; ipbrem < primPFBrem.size(); ipbrem++) {
0943 if (primPFBrem[ipbrem].indTrajPoint() == 99)
0944 continue;
0945 const reco::PFRecTrack& pfBremTrack = primPFBrem[ipbrem];
0946 double dist = pfBremTrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALShowerMax).isValid()
0947 ? LinkByRecHit::testTrackAndClusterByRecHit(pfBremTrack, clust, true)
0948 : -1.;
0949 if (dist > 0.) {
0950 if (nEcalDriven)
0951 iEnergy += clust.energy();
0952 else
0953 nEnergy += clust.energy();
0954 vecPFClusters.push_back(clust);
0955 }
0956 }
0957 }
0958 }
0959 }
0960 if (!vecPFClusters.empty()) {
0961 for (unsigned int pf = 0; pf < vecPFClusters.size(); pf++) {
0962 bool isCommon = ClusterClusterMapping::overlap(vecPFClusters[pf], *scRef);
0963 if (isCommon) {
0964 isSharingEnergy = true;
0965 }
0966 break;
0967 }
0968 }
0969 } else {
0970
0971
0972 bothGsfTrackerDriven = true;
0973 vector<PFCluster> nPFCluster;
0974 vector<PFCluster> iPFCluster;
0975
0976 nPFCluster.clear();
0977 iPFCluster.clear();
0978
0979 for (PFClusterCollection::const_iterator clus = theEClus.begin(); clus != theEClus.end(); clus++) {
0980 PFCluster clust = *clus;
0981 clust.calculatePositionREP();
0982
0983 float ndeta = fabs(nGsfPFRecTrack.gsfTrackRef()->eta() - clust.position().eta());
0984 float ndphi = fabs(nGsfPFRecTrack.gsfTrackRef()->phi() - clust.position().phi());
0985 if (ndphi > TMath::Pi())
0986 ndphi -= TMath::TwoPi();
0987
0988
0989 if (ndeta < 0.5 && fabs(ndphi) < 1.0) {
0990 double distGsf = nGsfPFRecTrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALShowerMax).isValid()
0991 ? LinkByRecHit::testTrackAndClusterByRecHit(nGsfPFRecTrack, clust)
0992 : -1.;
0993 if (distGsf > 0.) {
0994 nPFCluster.push_back(clust);
0995 nEnergy += clust.energy();
0996 } else {
0997 const vector<PFBrem>& primPFBrem = nGsfPFRecTrack.PFRecBrem();
0998 for (unsigned ipbrem = 0; ipbrem < primPFBrem.size(); ipbrem++) {
0999 if (primPFBrem[ipbrem].indTrajPoint() == 99)
1000 continue;
1001 const reco::PFRecTrack& pfBremTrack = primPFBrem[ipbrem];
1002 double dist = pfBremTrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALShowerMax).isValid()
1003 ? LinkByRecHit::testTrackAndClusterByRecHit(pfBremTrack, clust, true)
1004 : -1.;
1005 if (dist > 0.) {
1006 nPFCluster.push_back(clust);
1007 nEnergy += clust.energy();
1008 break;
1009 }
1010 }
1011 }
1012 }
1013
1014 float ideta = fabs(iGsfPFRecTrack.gsfTrackRef()->eta() - clust.position().eta());
1015 float idphi = fabs(iGsfPFRecTrack.gsfTrackRef()->phi() - clust.position().phi());
1016 if (idphi > TMath::Pi())
1017 idphi -= TMath::TwoPi();
1018
1019
1020 if (ideta < 0.5 && fabs(idphi) < 1.0) {
1021 double dist = iGsfPFRecTrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALShowerMax).isValid()
1022 ? LinkByRecHit::testTrackAndClusterByRecHit(iGsfPFRecTrack, clust)
1023 : -1.;
1024 if (dist > 0.) {
1025 iPFCluster.push_back(clust);
1026 iEnergy += clust.energy();
1027 } else {
1028 vector<PFBrem> primPFBrem = iGsfPFRecTrack.PFRecBrem();
1029 for (unsigned ipbrem = 0; ipbrem < primPFBrem.size(); ipbrem++) {
1030 if (primPFBrem[ipbrem].indTrajPoint() == 99)
1031 continue;
1032 const reco::PFRecTrack& pfBremTrack = primPFBrem[ipbrem];
1033 double dist = LinkByRecHit::testTrackAndClusterByRecHit(pfBremTrack, clust, true);
1034 if (dist > 0.) {
1035 iPFCluster.push_back(clust);
1036 iEnergy += clust.energy();
1037 break;
1038 }
1039 }
1040 }
1041 }
1042 }
1043
1044 if (!nPFCluster.empty() && !iPFCluster.empty()) {
1045 for (unsigned int npf = 0; npf < nPFCluster.size(); npf++) {
1046 for (unsigned int ipf = 0; ipf < iPFCluster.size(); ipf++) {
1047 bool isCommon = ClusterClusterMapping::overlap(nPFCluster[npf], iPFCluster[ipf]);
1048 if (isCommon) {
1049 isSharingEnergy = true;
1050 break;
1051 }
1052 }
1053 if (isSharingEnergy)
1054 break;
1055 }
1056 }
1057 }
1058
1059 return isSharingEnergy;
1060 }
1061 bool PFElecTkProducer::isInnerMost(const reco::GsfTrackRef& nGsfTrack,
1062 const reco::GsfTrackRef& iGsfTrack,
1063 bool& sameLayer) {
1064
1065
1066
1067 const reco::HitPattern& gsfHitPattern1 = nGsfTrack->hitPattern();
1068 const reco::HitPattern& gsfHitPattern2 = iGsfTrack->hitPattern();
1069
1070
1071 int gsfHitCounter1 = 0;
1072 for (auto const& hit : nGsfTrack->recHits()) {
1073 if (hit->isValid())
1074 break;
1075 gsfHitCounter1++;
1076 }
1077
1078 int gsfHitCounter2 = 0;
1079 for (auto const& hit : iGsfTrack->recHits()) {
1080 if (hit->isValid())
1081 break;
1082 gsfHitCounter2++;
1083 }
1084
1085 uint32_t gsfHit1 = gsfHitPattern1.getHitPattern(HitPattern::TRACK_HITS, gsfHitCounter1);
1086 uint32_t gsfHit2 = gsfHitPattern2.getHitPattern(HitPattern::TRACK_HITS, gsfHitCounter2);
1087
1088 if (gsfHitPattern1.getSubStructure(gsfHit1) != gsfHitPattern2.getSubStructure(gsfHit2)) {
1089 return (gsfHitPattern2.getSubStructure(gsfHit2) < gsfHitPattern1.getSubStructure(gsfHit1));
1090 } else if (gsfHitPattern1.getLayer(gsfHit1) != gsfHitPattern2.getLayer(gsfHit2)) {
1091 return (gsfHitPattern2.getLayer(gsfHit2) < gsfHitPattern1.getLayer(gsfHit1));
1092 } else {
1093 sameLayer = true;
1094 return false;
1095 }
1096 }
1097 bool PFElecTkProducer::isInnerMostWithLostHits(const reco::GsfTrackRef& nGsfTrack,
1098 const reco::GsfTrackRef& iGsfTrack,
1099 bool& sameLayer) {
1100
1101 unsigned int nLostHits = nGsfTrack->missingInnerHits();
1102 unsigned int iLostHits = iGsfTrack->missingInnerHits();
1103
1104 if (nLostHits != iLostHits) {
1105 return (nLostHits > iLostHits);
1106 } else {
1107 sameLayer = true;
1108 return false;
1109 }
1110 }
1111
1112
1113 void PFElecTkProducer::beginRun(const edm::Run& run, const EventSetup& iSetup) {
1114 auto const& magneticField = &iSetup.getData(magFieldToken_);
1115 auto const& tracker = &iSetup.getData(tkerGeomToken_);
1116
1117 mtsTransform_ = MultiTrajectoryStateTransform(tracker, magneticField);
1118
1119 pfTransformer_ = std::make_unique<PFTrackTransformer>(math::XYZVector(magneticField->inTesla(GlobalPoint(0, 0, 0))));
1120
1121 TransientTrackBuilder thebuilder = iSetup.getData(transientTrackToken_);
1122
1123 convBremFinder_ = std::make_unique<ConvBremPFTrackFinder>(thebuilder,
1124 mvaConvBremFinderIDBarrelLowPt_,
1125 mvaConvBremFinderIDBarrelHighPt_,
1126 mvaConvBremFinderIDEndcapsLowPt_,
1127 mvaConvBremFinderIDEndcapsHighPt_);
1128 }
1129
1130
1131 void PFElecTkProducer::endRun(const edm::Run& run, const EventSetup& iSetup) {
1132 pfTransformer_.reset();
1133 convBremFinder_.reset();
1134 }
1135
1136
1137 #include "FWCore/Framework/interface/MakerMacros.h"
1138 DEFINE_FWK_MODULE(PFElecTkProducer);